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Abstract 


Objectives :  Fugitive  dust  from  unconfined  sources  such  as  unpaved  roads  and  construction  areas 
can  impair  both  health  and  visibility.  In  this  report,  the  potential  benefits  of  vegetative 
windbreaks  in  reducing  vehicle-generated  fugitive  dust  are  studied  and  quantified  using  a 
combined  experimental  and  computational  approach.  The  overall  hypothesis  for  this  study  is  that 
maintaining  native  vegetation,  establishing  compatible  plantings  along  roads,  or  constructing 
windbreaks  could  be  a  useful  dust  mitigation  technique  on  military  training  ranges.  This  work 
focuses  on  PMio  dust  particles  with  an  aerodynamic  diameter  of  10  pm  or  less.  The  overall 
Objective  of  this  project  has  been  to  address  important  fugitive  dust  issues  for  the  Resource 
Conservation  and  Climate  Change  program  area  by  developing  and  validating  a  proof-of-concept 
computational  model  for  designing  windbreaks  for  dust  mitigation. 

Technical  Approach :  The  approach  toward  achieving  this  objective  has  been  to  utilize  laboratory 
and  field  experiment  data  to  develop  improved  models.  Model  development  was  accomplished 
by  integrating  particle  dispersion  and  atmospheric  turbulence  theory  with  wind  tunnel  and  field 
experiment  data.  The  computational  model  for  simulating  near-road  fugitive  dust  transport  was 
integrated  into  the  widely  used  Quick  Urban  and  Industrial  Complex  (QUIC)  Dispersion 
Modeling  System.  Prior  to  this  project,  QUIC  had  been  applied  to  near-source  pollutant  transport 
around  solid  obstructions  and  this  work  extends  its  applicability  to  incorporate  a  novel  sub¬ 
model  to  simulate  the  effects  of  both  thin  and  deep  vegetative  canopies  on  wind  fields  and 
particle  transport. 

Results'.  Results  from  this  project  include  improved  understanding  of  the  importance  of 
turbulence  in  enhancing  deposition  of  PMio  dust  particles  onto  vegetation.  The  results  show  that 
deposition  of  PMio  onto  surfaces  of  all  orientations  is  important.  The  results  indicate  that 
previous  estimates  of  vegetation  deposition  may  have  been  too  low. 

Benefits'.  The  benefits  from  the  project  include:  (i)  the  development  of  a  simple  predictive 
empirical  model  that  can  be  used  by  practitioners/end-users  to  design  vegetative  windbreaks  for 
reducing  near  source  fugitive  dust  particulate,  (ii)  the  development  of  a  validated  computation 
tool  (an  extension  to  the  QUIC  modeling  system)  for  studying  and  predicting  deposition  to  thin 
windbreaks  and  deep  vegetative  canopies. 
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1.  Objective 


This  project  was  designed  to  address  fugitive  dust  concerns  for  the  Resource  Conservation  and 
Climate  Change  program  area  of  Strategic  Environmental  Research  and  Development  Program 
(SERDP)  by  developing  computational  algorithms  and  input  data  for  predicting  fugitive  dust 
transport  and  emission  fluxes  within  a  few  hundred  meters  of  a  road  source.  The  objective  of 
this  project  was  to  develop  a  predictive  planning  windbreak  model  for  mitigating  near-road 
fugitive  dust.  This  was  accomplished  through  the  development  of  parameterizations  using 
laboratory  wind  tunnel  and  field  experiment  data.  These  parameterizations  were  used  to  achieve 
improved  prediction  of  transport,  dispersion  and  deposition  by  modifying  the  existing  Quick 
Urban  and  Industrial  Complex  (QUIC)  Dispersion  Modeling  System  (Singh  et  al.  2008;  Singh  et 
al.  2011).  QUIC  has  been  co-developed  by  Los  Alamos  National  Laboratory  and  University  of 
Utah.  A  real-time  running  prototype  of  the  modeling  system  was  demonstrated  in  the  field  at 
Hanford,  WA.  The  proposed  computational  model  for  simulation  of  fugitive  dust  mitigation  by 
windbreaks  addressed  two  Resource  Conservation  and  Climate  Change  areas:  (a)  prediction  of 
fugitive  dust  transport  and  emission  fluxes  and  (b)  integration  of  meteorological  conditions  and 
local  field  conditions  to  quantify  the  potential  transport  of  dust  plumes.  This  project  extended 
existing  computational  models  to  provide  a  proof-of-concept  tool  for  use  in  evaluating  the 
potential  dust  emissions  of  upcoming  training  activities  based  on  site-specific  terrain  and 
meteorological  data.  The  tool  can  also  be  utilized  in  mitigation  strategies. 

The  parameterizations  utilized  to  improve  QUIC  as  a  prediction  tool  included:  a  dry-deposition 
model  that  incorporated  the  enhancement  of  deposition  due  turbulence  (this  effect  had  been 
neglected  in  previous  models,  see  e.g.  Raupach  et  al.  2001)  and  models  for  the  influence  of 
windbreaks  on  mean  and  turbulent  wind  fields.  These  models  are  applicable  for  windbreaks  of 
varying  porosities  and  upstream  roughness  conditions. 

In  addition,  as  part  of  this  project,  we  used  our  improved  QUIC  simulations  and  field  data  to 
develop  a  simple  “empirical  equation”  predictive  model  that  can  be  used  by  practitioners/end- 
users  to  estimate  how  different  vegetation  configurations  can  reduce  near-source  emissions  by  a 
specified  fraction.  This  model  can  be  utilized  with  relatively  easy  to  obtain  infonnation  about  a 
specific  site. 

The  final  major  contribution  of  this  project  was  the  implementation  of  our  wind  tunnel  based 
parameterization  into  QUIC.  The  new  results  compare  well  with  the  experimental  data; 
reproducing  the  trends  and  the  magnitude  of  the  results  (Relative  errors  in  concentration 
predictions  for  the  two  field  experiments  were  12%  and  9.6%).  In  addition,  the  new  model 
strongly  indicates  that  the  reason  previous  dispersion  studies  required  such  unrealistically  small 
deposition  length  scales  for  vegetation  to  obtain  deposition  amounts  consistent  with  experimental 
data  was  the  failure  to  include  turbulence  in  the  deposition  parameterization. 
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2.  Background 

2.1  Description  of  the  Problem 

The  focus  of  this  project  has  been  on  trying  to  understand  the  interaction  between  vegetation  and 
particulate  matter  near  a  fugitive  dust  source  where  the  dust  plume  is  moving  horizontally 
through  vegetation.  We  have  concentrated  on  the  Environmental  Protection  Agency’s  criteria 
pollutant  PMio  (particulate  matter  with  a  diameter  of  10pm  or  smaller).  Near-road  fugitive  dust 
can  interact  with  thin  windbreaks  or  native  vegetation  canopies;  hence,  we  have  focused  our 
study  on  these  two  limiting  cases.  The  thickness  of  a  windbreak  is  generally  the  same  order  of 
magnitude  as  the  windbreak  height,  Hw,  or  less.  That  is,  windbreaks  are  relatively  “thin.”  Native 
vegetation  generally  extends  a  long  distance  downwind  of  a  road  source  and  consequently  has  a 
fetch  depth  that  is  very  large.  The  physics  of  the  flow  field  governing  these  two  situations  is 
quite  different.  The  contrasting  features  for  a  canopy  and  windbreak  are  shown  in  Figure  1.  Note 
that  the  coordinate  system  origin  is  located  along  the  vehicle’s  path  of  travel.  The  distance  from 
the  path  of  travel  to  the  windbreak  or  native  vegetation  is  taken  to  be  small  in  this  work.  The 
distance  between  the  road  and  the  start  of  the  vegetation  is  an  important  parameter  that  can 
substantially  impact  the  effectiveness  of  the  vegetation’s  ability  to  reduce  downwind 
concentrations.  The  further  from  the  road  the  vegetation,  the  more  time  the  dust  has  to  mix 
vertically,  effectively  increasing  the  “source  height”.  As  will  shown  in  the  results,  the  ratio  of  the 
source  height  to  the  height  of  the  vegetation  is  an  important  governing  parameter. 
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Figure  1:  Illustration  of  the  differences  between  windbreaks  and  canopies  along  with  the  coordinate 
system  used  for  simulations  and  field  studies  involving  unpaved  roads.  Windbreaks  are  very  thin  with  a 
fetch  depth  comparable  or  less  than  their  height,  while  canopies  have  an  infinitely  deep  fetch.  Both  are 
shown  in  the  context  of  the  current  project  downwind  of  an  unpaved  road. 

Referring  to  Figure  1,  assuming  that  the  road  is  long  and  that  end  effects  can  be  neglected, 
fugitive  dust  emissions  and  concentrations  are  independent  of  the  y-axis  because  of  the  road’s 
length.  Hence,  fugitive  dust  emissions  depend  only  upon  two  directions  (height  above  ground 
and  downwind  distance  from  the  road)  and  the  problem  can  be  assumed  to  be  quasi  two- 
dimensional. 

As  illustrated  in  Figure  2,  away  from  the  leading  edge,  the  flow  within  a  canopy  is  nearly  does 
not  substantially  change  with  increasing  downwind  distance.  Because  of  this,  the  flow  depends 
only  up  the  height  above  the  ground,  z.  In  contrast,  for  a  windbreak,  the  airflow  is  continually 
recovering  from  the  disturbance  caused  by  the  windbreak.  This  recovery  can  be  delayed  as  far 
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downwind  as  ~50  Hwb.  Because  of  this,  a  windbreak  profde  is  a  function  of  height,  z,  and 
downwind  distance,  x.  In  this  project,  we  have  used  several  simplifications  to  characterize 
important  physical  vegetative  characteristics  that  effect  both  momentum  and  particle  transport. 
For  example,  the  behavior  of  the  flows  for  both  windbreaks  and  canopies  depend  on  the  amount 
of  vegetative  surface  area  present.  For  canopies,  this  can  be  expressed  by  Leaf  Area  Index  (LAI). 
LAI  is  defined  as  the  vegetative  surface  area  per  unit  of  ground  covered  by  the  canopy.  In  this 
work  we  assume  that  vegetative  surface  area  is  distributed  uniformly  through  the  canopy  (both 
vertical  and  horizontally).  The  density  of  the  vegetative  or  y  is  obtained  by  dividing  the 
vegetative  surface  area  by  the  volume  of  the  canopy.  It  is  the  surface  area  per  unit  volume.  For 
windbreaks  optical  porosity,  P,  is  the  measureable  parameter  indicating  thickness,  y  can  be 
obtained  from  P  by  using  the  relationships  from  the  literature  (e.g.  Raupach  et  al.  2001). 


Windbreak 


u{  z) 


Figure  2:  Velocity  profiles  further  illustration  of  the  differences  between  windbreaks  (top  panel)  and 
canopies  (lower  panel).  The  mean  velocity  profile  downstream  of  a  windbreak  is  continuously  changing 
and  slowly  approaches  the  upstream  profile,  while  a  mean  canopy  profile  is  relatively  independent  of 
position  within  the  canopy. 
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2.1  Windbreak  Literature 

There  have  been  many  studies  focusing  on  different  aspects  of  windbreaks  (or  shelterbelts)  and 
their  applications  in  the  literature.  Windbreak  Technology  from  the  Proceedings  of  an 
International  Symposium  on  Windbreak  Technology  (Brandle,  1988)  provides  an  excellent 
overview  of  the  many  uses  of  windbreaks  including  their  role  in  mitigating  wind  erosion,  snow 
movement  and  agricultural  spray  drift;  their  role  changing  crop  yield,  as  well  as  surface  energy 
balances  and  evaporation  rates.  Many  field  campaigns,  wind  tunnel  experiments,  analytical 
analysis,  and  numerical  modeling  activities  have  been  undertaken  in  this  field.  Field  studies  of 
windbreaks  have  included  the  work  of:  Bradley  and  Mulheam  (1983),  Nord  (1991),  and  Wilson 
(2004).  Wind  tunnel  studies  have  included:  Plate  (1971),  Judd  et  al.  (1996),  Perrera  (1981),  and 
Guan  et  al.  (2003).  Numerical  modeling  studies  have  included  using  various  fonns  of  the 
Reynolds  Averaged  Navier-Stokes  equations  (Wang  et  al.,  2001;  Wilson,  1985;  Bourdin  and 
Wilson,  2008;  and  Santiago,  2007).  In  addition,  high-resolution  Large-Eddy  Simulations  of 
windbreaks  were  conducted  by  Patton  et  al.  (1998).  Analytical  analysis  include:  Counihan  et  al. 
(1974)  and  Wilson  (1990).  A  few  of  these  studies  develop  empirically  based  expressions  for  the 
mean  wind  and  turbulence  (some  of  which  have  been  utilized  in  the  present  work).  These  studies 
include:  Perrera  (1981),  Bradley  and  Mulheam  (1983),  Counihan  et  al.  (1974),  and  Guan  (2003). 
In  the  work  conducted  for  this  SERDP  project,  we  have  drawn  extensively  from  this  literature. 

Despite  the  common  need  of  the  current  work  and  other  windbreak  studies  to  either  model  or 
characterize  the  mean  and  turbulence  wind  fields,  there  are  substantial  differences  with  respect  to 
the  particles  studied.  Tabler  (1991  and  1992)  examined  the  mitigation  of  transport  of  snow 
across  highways.  Other  windbreak  studies  have  been  concerned  with  the  transport  of  sand 
(Musick  and  Gillette,  1990),  which  typically  has  a  size  range  between  60-2000  pm.  Another 
important  large  particle  application  is  spray  drift  mitigation  of  pesticides/insecticides,  which 
have  droplet  diameters  that  are  several  hundred  microns  (Wood  et  al.,  2001).  Such  particles  are 
much  larger  than  PMio.  Because  of  the  larger  particle  sizes  appearing  in  this  part  of  the 
windbreak  literature,  the  relevant  removal  mechanisms  are  different  from  the  current  work.  PMio 
tends  to  be  impaction  dominated  rather  gravitational  settling  dominated,  which  is  usually  the 
case  for  the  large  particle  problems.  Gravitational  settling  is  well  understood  and  accurately 
described  by  parameterizations  in  the  literature  (e.g.  Seinfeld  and  Pandis,  1998).  Impaction  is 
much  less  understood  and  a  significant  aspect  of  the  current  work  is  to  improve  impaction 
modeling  for  PMi0.  Since  the  primary  deposition  mechanisms  associated  with  large  particle  and 
small  particle  deposition  are  fundamentally  different  (gravitational  settling  versus  impaction),  the 
QUIC  simulations  for  the  transport  of  PMio  are  not  comparable  to  the  transport  of  larger  particles 
such  as  sand.  There  is,  however,  a  need  to  tie  QUIC  and  the  findings  from  this  work  to  these 
other  large  particle  problems  since  many  applications  may  have  a  wide  range  of  particle  sizes. 
This  is  discussed  further  in  the  recommendations  for  future  work. 

In  this  work,  we  have  attempted  to  experimentally  characterize  the  critical  factors  influencing 
particle  transport  and  deposition  in  vegetation  so  that  modeling  of  these  processes  can  be 
improved.  In  the  following  sections  we  first  present  our  technical  approach  (methodology  and 
theory)  followed  by  important  results  from  our  research.  The  technical  approach  and  results  are 
organized  by  basic  tasks  that  needed  to  be  completed.  Hence,  we  begin  with  a  description  of  the 
wind-tunnel  studies  that  were  used  to  improve  our  understanding  of  turbulence  on  particle 
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deposition  and  to  develop  a  new  deposition  model.  We  then  transition  to  numerical  model 
development  for  the  extended  QUIC  model  and  the  theoretical  foundations  of  the  sub-models. 
Next,  a  discussion  of  the  new  field  experiments  that  were  conducted  during  this  project  is 
presented.  Finally,  a  description  of  the  new  simple  “empirical  equation”  predictive  model  is 
presented. 

3.  Materials  and  Methods 

The  project  can  be  broken  up  in  to  three  areas:  laboratory  measurements,  field  experiments,  and 
model  development  and  prototype  demonstration.  In  the  following  sections,  we  begin  by 
describing  the  laboratory  measurements,  then  the  modeling  approach,  followed  by  the  field 
experiments. 

3.1.  Laboratory  Wind  Tunnel  Measurements 

3.1.1.  Deposition  Experiments 

Wind  tunnel  experiments  were  perfonned  to  test  our  hypothesis  that  enhanced  turbulent  mixing 
from  the  interaction  of  a  vegetative  canopy  with  vehicle-generated  dust  will  increase  total 
deposition  onto  vegetative  surfaces.  More  specifically,  several  existing  models  are  based  on 
laminar  flow  impaction  theory  (e.g.  Raupach  et  ah,  2001;  Aylor  and  Flesch,  2001)  for  deposition 
to  vegetative  canopies  and  windbreaks  that  do  not  account  for  the  details  of  the  turbulence  eddy 
sizes  and  strength  of  the  local  fluctuations.  Previous  field  studies  have  reported  and  observed 
increased  deposition  in  a  vegetative  canopy  (Veranth  et  al.  2003).  Deposition  measurements 
were  accomplished  using  a  fluorimetry  measurement  technique  that  measures  known  amounts  of 
fluorescein  (C2oHi0Na205,  Sigma-Aldrich)  aerosol  particles  deposited  onto  substrates.  Three 
different  sizes  (0.25  cm2,  1.00  cm2,  and  1.44  cm2)  of  hard  plastic  polypropylene  square  substrates 
(1  mm  thick)  were  attached  on  a  six-axis  deposition  frame.  Small  holes  were  drilled  along  the 
centerline  of  the  frame  and  24-gauge  wire  was  woven  into  the  frame  to  create  a  directional 
matrix  for  the  substrates  to  be  placed.  Substrates  were  held  onto  wire  using  a  pressure-sensitive 
adhesive  (Blu-Tack,  Bostick  Inc.).  Figure  3  shows  a  schematic  of  the  experimental  apparatus 
used. 
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Test  Section 


Figure  3:  Schematic  of  the  wind  tunnel  deposition  experimental  apparatus  (not  to  scale). 


To  test  our  hypothesis  of  deposition  enhancement  from  turbulent  motions,  a  grid  designed  to 
generate  isotropic  turbulence  was  constructed  with  a  solidity  or  open  area  ratio  of  0.5  and  placed 
just  downstream  of  the  converging  section  of  the  wind  tunnel  upstream  of  the  deposition  surfaces. 
Grid-generated  turbulence  was  selected  for  these  experiments  because  it  may  be  used  as  a  close 
approximation  to  idealized  homogenous  isotropic  turbulence.  Since  this  type  of  flow  has  been 
studied  extensively  both  experimentally  and  theoretically  (see  Pope  2000  for  a  review),  it 
provides  an  ideal  experimental  framework  for  relating  various  scales  of  turbulence  to  deposition 
characteristics. 

The  aerosol  generator  used  was  an  Ultra  Sonic  Humidifier  (Model  V5100NS,  Kaz  Inc.).  It  was 
positioned  on  a  stand  beneath  the  wind  tunnel  and  the  aerosol  was  piped  up  to  the  centerline  of 
the  wind  tunnel  using  polyethylene  tubing.  The  aerosol  injected  into  the  wind  tunnel  was  a 
mixture  of  80%  distilled  water,  20%  glycerol,  and  0.5  mg  of  soluble  fluorescein  per  ml  of 
solution.  A  laser-based  aerosol  spectrometer  (Model  1.109  Grimm  Technologies  Inc.)  was  used 
to  analyze  the  size  range  and  concentration  of  aerosol  being  measured  at  the  points  of  deposition 
inside  the  wind  tunnel.  The  aerosol  diameter  (after  initial  water  evaporation)  was  relatively 
evenly  distributed  over  the  size  range  between  1-5  pm. 

Regardless  of  substrate  size,  location,  or  mean  wind  velocity,  the  aerosol  was  injected  into  the 
test  section  for  one  hour  to  allow  for  accurate  mass  deposition  measurements.  Each  experiment 
was  repeated  three  times  as  the  substrate  locations  were  rotated  clockwise  after  each  experiment. 
The  procedure  was  used  to  average  out  any  variations  or  errors  that  might  exist  in  the 
concentration  amounts  at  individual  substrate  locations  and  to  carry  out  each  experiment  on  all 
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six  axes  simultaneously.  Figure  4  shows  a  schematic  for  how  the  substrates  were  placed  onto  the 
deposition  frame.  Substrate 


Side  View 


Front  View 


Figure  4:  Side  view  (left)  and  front  view  (right)  of  the  deposition  frame  showing  plastic  substrates 
attached  to  24-gauge  wire.  Note  in  the  side  view,  the  flow  is  from  left  to  right  and  in  the  front  view  flow 
is  directed  into  the  page. 


locations  were  designated  in  accordance  to  the  right-hand  rule  facing  downstream  of  the  wind 
tunnel.  Since  6  substrates  were  used  (±x,  ±y,  ±z)  for  any  given  experiment,  the  substrates  were 
sandwiched  together  using  Blu-Tack,  while  the  24-gauge  wire  was  used  as  the  structure  to  keep 
the  substrates  in  a  fixed  location  without  substantially  modifying  the  flow  field.  The  amount  of 
fluorescein  deposited  onto  each  substrate  was  measured  using  a  fluorimeter,  which  enabled 
quantification  of  the  deposition  fraction  ( DF ).  Here,  DF  is  defined  as  the  ratio  of  the  mass  of 
fluorescein  deposited  on  each  substrate  to  the  total  amount  of  fluorescein  available  in  the  volume 
of  air  intercepted  by  the  substrate  cross  section  over  one  hour.  In  order  to  obtain  different 
turbulence  intensities  and  length  scales,  the  deposition  frame  was  placed  at  1.80  and  2.13  m 
downstream  of  the  grid. 


After  each  deposition  experiment,  the  substrates  were  carefully  removed  from  the  wires  using 
surgical  clamp  forceps,  and  then  washed  in  a  test  tube  with  2  ml  of  distilled  water.  Each  substrate 
was  individually  washed  in  distilled  water  and  analyzed  in  a  fluorimeter.  The  fluorimeter  outputs 
a  dimensionless  reading  of  fluorescence  intensity.  Therefore,  in  order  to  quantify  the  amount  of 
fluorescein  that  was  deposited  on  each  substrate  a  dilution  series  was  used  to  develop  a 
calibration  curve.  A  linear  regression  analysis  was  used  to  correlate  concentration  amounts  of 
fluorescein  deposited  on  individual  substrates  to  the  output  fluorescence  measured  by  the 
fluorimeter.  A  good  linear  fit  (R2  =0.998)  was  obtained  over  the  range  of  experimental  data.  The 
linear  fit  was  then  used  to  calculate  the  mass  of  fluorescein  that  was  deposited.  The  deposition 
fraction  was  determined  as  follows: 
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where,  M  is  the  ratio  of  fluorescein  mass  deposited  to  the  total  amount  of  available  fluorescein  in 
the  airstream  over  a  one  hours  sampling  period,  U0  is  the  mean  velocity,  A  is  the  area  of  the 
substrate  used,  t  is  the  run  time  of  the  experiment,  and  C  is  the  concentration  value  measured  at 
the  location  of  deposition.  The  concentration  was  determined  from  optical  measurements  using  a 
GRIMM  1.109  particle  spectrometer  instrument.  DF  was  calculated  for  all  six  axes. 

3.1.2.  Vegetation  Deposition  Experiments 

Wind  tunnel  deposition  experiments  were  also  conducted  on  various  vegetative  elements 
including  artificial  broad  leaf,  branch,  grass,  and  metal  welding  rod  elements  (see  Table  1  for  a 
summary  of  experiments).  The  experimental  procedure  was  very  similar  to  the  above  artificial 
substrate  deposition  experiments  using  the  same  wind  tunnel,  aerosol  generator,  fluorescein 
mixture,  and  fluorimeter.  However,  unlike  the  artificial  substrate  experiments,  the  deposition 
experiments  were  only  conducted  at  1.80  m  downstream  and  only  with  the  turbulence  grid  in 
place.  The  purpose  of  these  experiments  was  to  verify  that  the  turbulence  enhanced  deposition 
parameterizations  found  for  the  artificial  substrates  were  valid  for  more  realistic  vegetative 
elements.  The  primary  difference  between  the  vegetation  experiments  and  the  artificial  substrate 
experiments  is  that  the  resulting  deposition  fraction  is  an  integral  quantity  over  the  entire  area  of 
the  vegetative  element.  The  fluorescein  mass  deposited  was  measured  for  the  entire  element 
rather  than  isolating  for  each  surface  and  orientation.  It  is  believed  that  this  better  matches  what 
would  be  required  for  a  real-world  vegetative  deposition  efficiency  model. 

The  area  of  the  irregular  leaf,  grass,  and  twig  elements  was  required  to  calculate  the  deposition 
efficiency  (Eq.  1).  In  the  case  of  the  broad  leaf  and  grass  samples,  direct  measurement  of  the  area 
was  impractical.  A  pixel  counting  technique  was  used  to  measure  the  frontal  area.  Each  of  the 
samples  was  photographed  with  a  Canon  XSi  DSLR  camera  on  grid  paper  and  the  photo  was 
loaded  into  Adobe  Photoshop  CS5.5.  Using  the  “Magnetic  Lasso”  tool,  the  perimeter  was  traced 
via  contrast  detection  methods.  The  number  of  pixels  in  the  leaf  area  and  known  area  of  grid 
cells  on  the  paper  were  recorded.  The  area  of  the  leaf  was  then  calculated  by  Eq.  2. 


Leaf  Area  = 


Control  Area 


x  No.  of  Leaf  Pixels  . 


No. Pixels  in  Control  Area 


(2) 
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Substrate 

Orientation 

Average  Max 
Width  or 
Diameter  (mm) 

Average  Area  (mm2) 

Square 

Polypropylene 

±  X,  Y,  Z 

5,  10,  12 

25, 100, 144 

Broad  Leaf 

0°,  +45°,  - 
45°,  Edge -On 

26.3 

1279.6 

Onion  Grass 

0° 

1.0 

47.0 

Twig 

0° 

2.5 

192.0 

Welding  Rods 

0° 

0.89,  1.57,2.39 

44.45,  78.74,  119.38 

Table  1:  Summary  table  of  the  various  deposition  experiments  run  in  the  wind  tunnel.  Note  that  a  0° 
orientation  indicates  that  vegetative  surface  of  greatest  area  was  oriented  normal  to  the  incoming  wind. 


Broad  leaf  and  grass  area  were  calculated  with  this  method  and  it  was  found  to  be  very 
repeatable  with  <1%  error.  The  area  of  the  twigs  was  perfonned  with  a  statistical  sampling 
technique  where  the  diameter  was  measured  at  6  positions  along  the  length  and  the  area  was 
calculated  using  the  average  diameter  measurement  and  the  length.  The  projected  frontal  area  for 
the  welding  rods  was  calculated  by  simply  measuring  the  diameter  and  multiplying  by  the  length. 


In  addition  to  taking  photographs  to  calculate  the  area  of  the  leaves  and  grass  samples,  pictures 
were  taken  of  the  samples  before,  during,  and  after  the  wind  tunnel  deposition  experiments. 
Photographs  were  taken  of  the  leaf  (or  grass)  in  position  without  any  air  flow,  just  after  initial 
startup,  near  the  end  of  the  test,  and  after  the  tunnel  had  completely  stopped.  The  purpose  of 
these  pictures  was  to  quantify  the  streamlining  of  the  leaf  in  the  wind  and  any  permanent 
relaxation  after  the  test  was  finished. 


Due  to  the  various  shapes  and  deposition  regimes  tested  in  the  wind  tunnel,  the  run  time  required 
for  each  element  type  to  achieve  sufficient  fluorescence  signal  without  over  saturating  the 
surface  varied  according  to  element  width  and  the  tunnel  wind  speed.  Also,  the  volume  of  the 
washing  fluid  varied  from  1-10  mL  based  on  the  physical  size  of  the  element. 


3.1.3.  Hot-wire  Turbulence  Quantification  Experiments 

A  single-sensor  straight-wire  hot-wire  probe  (55P16  Probe,  Dantec  Dynamic  Inc.)  was  placed  in 
the  wind  tunnel  to  gather  turbulence  data  to  correlate  deposition  enhancement  to  the  proper 
scales  of  turbulence  in  the  wind  tunnel.  Traditionally,  the  Taylor-scale  Reynolds  (R>  )  number  is 
used  to  characterize  grid  turbulence  (Pope  2000),  therefore,  hot-wire  probes  were  used  to 
accurately  determine  the  streamwise  r.m.s.  of  the  streamwise  velocity  (u’)  and  the  Taylor  micro¬ 
scale  (kg): 


(3) 


We  hypothesize  that  the  Taylor  micro-scale  (kg)  is  the  turbulence  scale  most  responsible  for 
enhanced  deposition  because  this  length  scale  is  shown  to  have  the  smallest  dynamic  eddies 
(Tennekes  and  Lumley  1972).  At  these  smaller  scales,  we  expect  the  eddy-sizes  of  the  micro¬ 
scale  to  interact  with  the  aerosol  particles  and  contribute  in  actively  depositing  these  particles 
onto  the  substrate  surfaces.  These  eddies  are  still  capable  of  transferring  energy  to  much  smaller 
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dissipative  scales.  Therefore,  we  believe  Lg  should  play  an  important  role  in  characterizing 
deposition  onto  known  surfaces.  This  is  shown  in  Figure  5. 


Figure  5:  Schematic  showing  the  interaction  of  turbulent  eddies  with  a  size  on  the  order  of  the  Taylor 
micro-scale  with  the  particles  that  deposit  onto  the  surface  of  the  substrates. 


The  hot-wire  probe  was  calibrated  using  the  University  of  Utah  portable  hot-wire  calibration  test 
facility  before  each  experiment.  Data  were  sampled  at  25kHz  using  an  AN- 1003  Anemometry 
system  (AA  Lab  Systems  Ltd),  BNC  2110  I/O  ADC  (National  Instruments),  NI-6122  DAQ 
(National  Instruments),  and  Labview  10.  A  differential  pressure  transducer  (10  Torr  MKS 
Baratron)  and  hot-wire  probe  were  colocated  at  the  exit  of  the  converging  calibration  facility.  A 
third-order  polynomial  fit  was  applied  in  order  to  develop  a  calibration  curve  and  equation  for 
the  wind  tunnel  measurements  (Bruun  1995;  Tennekes  and  Lumley  1972).  Hot-wire  calibration 
was  perfonned  both  before  and  after  the  data  were  acquired  in  the  wind  tunnel  and  the  average  of 
the  two  calibration  equations  was  used  in  the  final  data  analysis.  The  calibration  is  shown  below 
in  Figure  6.  Data  were  then  taken  at  four  separate  locations  in  the  test  section  of  the  wind  tunnel. 
These  locations  were  identical  to  those  of  the  deposition  measurements. 


Figure  6:  Hot-wire  probe  velocity  (U)  -  voltage  (E)  calibration  curve.  Calibration  was  performed  both 
before  (subscript  i)  and  after  (subscript  f)  data  collection  in  the  wind  tunnel.  The  average  curve  (subscript 
ave)  was  used  to  process  the  data.  A,  B  and  n  represent  best-fit  coefficients. 
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3.2.  QUIC  Model  Development 

An  important  aspect  of  this  project  has  been  the  development  of  new  models  for  predicting 
deposition  onto  vegetative  windbreaks.  In  particular,  added  the  new  models  to  an  existing 
software  system  maintained  at  Los  Alamos  National  Laboratory  called  the  QUIC  (Quick  Urban 
and  Industrial  Complex)  Dispersion  Modeling  System-  (Pardyjak  and  Brown  2001;  Williams  et 
al.  2004;  Singh  et  al.  2008;  Singh  et  al.  2011).  QUIC  has  been  co-developed  by  researchers  at  the 
University  of  Utah  and  Los  Alamos  National  Laboratory  over  the  past  decade  and  computes  3D 
high-resolution  flow  fields  and  particle  dispersion  around  building  in  urban  areas.  QUIC  is 
composed  of  a  wind  solver  (QUIC-URB),  which  computes  wind  fields  based  using  initialized 
wind  fields  based  on  simple  experimental  parameterizations  and  then  forces  the  flow  field  to  be 
mass  conservative  (Singh  et  al.  2008).  For  dispersion  modeling  the  system  uses  QUIC-PLUME, 
a  random  walk  or  Lagrangian  particle  dispersion  model  (Williams  et  al.  2004).  QUIC  also  has  an 
easy  to  use  Graphic  User  Interface  (GUI)  that  allows  for  building  the  simulation  domain  and 
visualizing  results.  Prior  to  this  project,  QUIC  did  have  not  vegetative  windbreak  capabilities, 
however,  because  of  its  ability  to  simulate  urban  and  other  complex  flow  and  dispersion 
problems  and  its  free  availability,  it  represented  an  excellent  base  to  build  upon.  Figure  7  below 
shows  an  illustrative  example  of  the  type  of  results  that  are  easily  obtained  using  the  QUIC 
model. 


Figure  7:  Illustration  of  an  urban  flow  and  dispersion  simulation  using  QUIC.  QUIC  is  the  based  model 
that  was  built  upon  in  this  project. 
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3.3.1  Improved  Mean  Wind  and  Turbulence  Modeling  for  Windbreaks 

Windbreaks,  such  as  a  single  or  few  rows  of  trees  or  fences  are  common  roughness  features 
found  alongside  roads  (Judd  et  al.  1996).  A  goal  of  the  current  SERDP  project  is  to  successfully 
model  horizontal  advection  and  turbulent  diffusion,  along  with  deposition  within  single  and 
arrays  of  windbreaks  within  the  QUIC  platfonn. 

A  qualitative  description  of  flow  resulting  from  the  placement  of  windbreaks  is  shown  in  Figure 
8  and  is  documented  in  Judd  et  al.  (1996).  Starting  on  the  left  of  the  figure,  one  finds  the  typical 
upwind  approach  profile  found  in  the  atmospheric  surface  layer.  The  letter  A  indicates  this.  As 
the  flow  approaches  the  windbreak,  the  flow  slows  down  in  the  streamwise,  or  x  direction,  and  as 
a  consequence  of  mass  conservation,  the  vertical  velocity  component  increases.  On  the  leeward 
side  of  the  windbreak  the  flow  is  divided  into  two  regions  dependent  upon  height:  for  z>Hwb 
there  is  an  accelerated  displaced  flow,  denoted  by  the  letter  C.  For  z<Hwb  flow  velocities  and 
turbulence  intensities  are  lower  than  in  the  approach  profile.  The  flow  for  z<Hwb  is  divided  into 
the  bleed  flow  region  directly  leeward  of  the  windbreak,  and  the  quiet  zone  located  downwind  of 
the  bleed  flow  region.  These  zones  are  denoted  by  the  letters  B  and  D  respectively.  In  the  red 
cone  shaped  region  (shear  zone  denoted  by  letter  E),  the  displaced  flow  interacts  with  the  quiet 
zone  in  an  area  containing  relatively  high  shear.  As  the  shear  zone  entrains  fluid,  it  decreases  in 
intensity  and  eventually  contacts  the  ground  and  develops  into  a  re-equilibrium  zone  denoted  by 
the  letter  F.  In  the  re-equilibrium  zone  the  profile  slowly  develops  back  into  the  approach  profile. 


Figure  8:  Flow  zones  and  velocity  profiles  downwind  of  a  windbreak  as  defined  in  Judd  et  al.  (1996) 
together  with  the  road  orientation  utilized  in  the  current  simulations  and  field  studies. 

Mean  Wind 

This  section  details  the  methods  utilized  to  implement  the  fast  response  windbreak  model  within 
the  constructs  of  the  QUIC  dispersion  modeling  tool  platfonn.  To  summarize,  the  QUIC 
platform  utilizes  an  empirically  based  initial  wind  field.  This  initial  wind  field  will,  in  general, 
not  be  mass  conservative.  QUIC  utilizes  a  basic  variational  analysis  procedure  to  modify  the 
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initial  wind  field  into  a  mass  conservative  wind  field.  The  majority  of  the  research  effort  in  the 
current  project  is  to  define  this  empirically  based  initial  wind  field  (Singh  et  al.  2008). 

The  basic  modeling  procedure  used  here  utilizes  existing  high-resolution  wind  tunnel  data  sets. 
Judd  et  al.  (1996)  and  Guan  et  al.  (2003)  performed  wind  tunnel  studies  of  windbreaks  with 
varying  optical  porosities  and  proposed  some  simple  models  that  we  are  utilizing  for  the 
initial  non-mass  conservative  velocity  field  in  QUIC.  The  basic  methodology  is  illustrated  in 
Figure  9.  The  physical  region  upwind  and  downwind  of  an  infinitely  thin  windbreak  is  sub¬ 
divided.  A  different  parameterization  is  utilized  in  each  region  for  the  initial  non-  mass 
conserved  velocity  field.  For  Region  3,  shown  in  red,  we  adopt  the  similarity  solution 
parameterization  proposed  by  Perera  (1981).  This  solution  is  a  function  of  the  optical  porosity, 
windbreak  height,  downwind  distance,  height,  upstream  surface  roughness  and  displacement 
height.  It  is  a  similarity  solution  that  is  not  applicable  directly  downwind  of  the  windbreak.  It  is 
applicable  at  distances  at  least  7.5  times  the  windbreak  height  (7.5  Hwb )  downwind. 


Region  1 


UpwindProfile  (a„=l) 


Upwind  of 

< -  Wind  Break 


Region  2 


Tanh  Profile  (cu=tanh  (f(x),  z/H^bJ) 


Tanh  profiles 

Less  than  7.5  Wind 

_ | _  Break  Heights  down 

Wind 


A:  Approach  profile 


B 


C:  Displaced  profile 


E:  Mixing  zone  (Red  Zone) 


B:  Bleed  flow  D:  Quiet  zone  F:  Re-equilibrium  zone 

Figure  9:  Illustration  of  division  of  windbreak  flow  geometry  into  three  regions  that  extend  from  the 
ground  to  an  infinite  height  vertically.  The  border  between  Region  1  and  2  (the  blue  and  white  regions)  is 
the  location  of  the  windbreak.  The  border  of  Region  2  and  Region  3  (the  white  and  red  regions)  is  a  line 
extending  from  the  ground  vertically  to  the  top  of  the  simulation  domain  located  1.5Hwb  downwind  of  the 
windbreak. 


In  Region  1  (shown  in  blue),  upwind  of  the  windbreak,  changes  in  the  upstream  profile  are 
minimal.  However,  as  previously  noted  there  is  a  deceleration  in  the  profile  as  it  approaches  the 
windbreak,  as  well  as  acceleration  over  the  top  of  the  windbreak.  It  is  assumed  that  the 
variational  procedure  in  QUIC  will  be  sufficient  to  model  these  components.  Hence,  all 
locations  upwind  showed  in  blue  use  a  simple  boundary  layer  profile  for  the  initial  non-mass 
conserved  field.  For  Region  2,  shown  in  white,  the  aerodynamic  porosity  as  defined  by  Guan  et 
al.  (2003)  is  helpful.  It  is  defined  in  Eq.  4,  as: 
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(4) 


a(x,z) 


u(x,z) 

u„(x,z) 


Where  a(x,z )  is  the  aerodynamic  porosity,  ua(z )  is  the  upstream  velocity  at  height  z,  and  u(x,z)  is 
the  local  velocity,  a  is  1  when  the  local  velocity  and  the  upstream  velocity  are  identical.  Any 
reduction  in  momentum  will  decrease  a. 


In  Region  2,  we  propose  to  model  a  and  then  multiply  a  by  the  upstream  velocity  at  the  given 
height  to  determine  the  velocity,  rather  than  directly  estimating  the  velocity.  The  reason  for  this 
indirect  method  is  that  a  is  generally  easier  to  model  than  the  wind  velocity.  To  illustrate  this 
point,  take  for  example,  the  upwind  profile.  It  can  be  many  different  relatively  complex  curves. 
However,  the  upwind  aerodynamic  porosity  is  always  a  very  simple  constant  function  equal  to  1 . 

In  Region  2,  a  one-dimensional  vertical  function  is  implemented.  The  inputs  for  this  vertical 
function  vary  with  increasing  distance  downwind  from  the  windbreak.  The  behavior  of  the 
function  needed  in  the  direct  downwind  zone  is  as  follows:  at  heights  much  greater  than  the 
windbreak,  for  example  2  to  3  times  Hwb  a  is  1.  At  heights  just  behind  the  fence,  for  z<H„b,  a  is  a 
function  of  the  windbreak’s  optical  porosity.  Also,  there  is  a  similarity  of  the  flow  directly 
downwind  of  the  windbreak  to  that  of  a  classical  mixed  layer.  These  properties  suggest  that  a 
hyperbolic  tangent  function,  tank,  would  be  appropriate  to  model  a  in  Region  2  just  downwind 
of  the  windbreak.  The  tank  function  is  a  well-behaved  function  with  limits  of  1  and  - 1  that  are 
approached  as  the  function  inputs  approach  oo  and  -oo  respectively.  High  above  the  windbreak  the 
function  naturally  approaches  the  limit  for  a,  of  1 .  Near  the  ground,  behind  the  windbreak,  the 
function  may  be  modified  to  approach  the  aerodynamic  porosity  of  the  windbreak.  Visually  this 
is  shown  in  Figure  9.  With  the  mathematical  function  defined,  two  sub  models  are  needed  to 
complete  the  zone  and  are  described  in  following  section. 

Windbreak  Bleed  Flow  Model 

Guan  et  al.  (2003)  propose  a  model  for  the  flow  through  a  windbreak,  the  so-called  “bleed  flow.” 
A  bleed  flow  model  relates  a  windbreak’s  optical  porosity,  P,  to  the  value  of  a  directly  behind  the 
fence.  Guan  et  al.  (2003)  propose  two  different  models:  1)  for  thin  windbreaks  such  as  fences 
and  2)  for  relatively  thick  fences  such  as  rows  of  trees.  The  models  are  presented  respectively  as: 

a  =  p  (thin  windbreaks)  (5) 

a  =  ffA  (thick  windbreaks)  (6) 

The  value  of  a  then  replaces  the  lower  limit  of  -1  for  tanh  as  the  ground  is  approached. 
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Figure  10:  A  detailed  view  of  Region  2.  Illustrations  of  the  use  of  the  tank  function  in  vertical  velocity 
profdes  on  the  leeward  side  of  the  windbreak.  Note  the  origin  of  the  tank  function,  the  dotted  horizontal 
line  originating  at  the  top  of  the  windbreak.  Note  also  the  top  and  bottom  boundaries  of  the  shear  zone. 
The  difference  between  these  two  lines  is  the  shear  zone  depth  8. 

Windbreak  Shear  Wake  Model 

An  additional  step  in  defining  the  zone  adjacent  and  downwind  of  the  windbreak  is  the  rate  of 
spread  of  a  zone  of  relatively  high  shear  just  downwind  of  the  top  of  the  windbreak.  These  are 
shown  in  Figure  10  by  the  blue  dotted  lines  than  run  diagonally  downwind  of  the  top  of  the 
windbreak  and  are  labeled  “top  boundary  of  shear  zone”  and  “bottom  boundary  of  shear  zone.” 
Judd  et  al.  (1996)  present  the  following  model  to  estimate  the  rate  of  spread  of  this  high  shear 
zone  in  the  stream  wise  direction: 
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Where,  UHwb  is  the  upwind  velocity  at  the  windbreak  height,  aw  is  the  variation  of  the  vertical 
velocity  component,  w,  A u  is  uHwb*(l-a),  and  uave  is  UHwb*(l  +a)/2.  This  model  defines  the  rate  of 
spread  of  the  shear  zone.  The  thickness  of  the  shear  zone  is  0  at  the  top  of  the  windbreak.  It  then 
increases  as  the  Judd  et  al.  (1996)  model  dictates.  To  relate  the  rate  of  spread  to  the  tank 
function,  a  definition  utilized  by  researchers  of  classical  mixing  layers  is  needed:  the  thickness 
of  this  shear  zone  accounts  for  90%  of  the  difference  between  the  two  velocity  streams.  This 
90%  criterion  corresponds  to  inputs  of  1.5  and  -1.5  in  the  tank  function. 

Lastly,  the  origin  of  the  tank  function  must  be  detennined.  This  is  shown  as  a  horizontal  dotted 
line  originating  at  the  top  of  the  fence  and  continuing  on  the  leeward  side  in  Figure  10.  Judd  et  al. 
(1996)  note  that  the  inflection  point  of  the  velocity  profile,  the  origin  of  the  tank  system,  is 
nearly  located  on  one  streamline  that  passes  through  the  pointed  located  at  the  windbreak  height 
in  the  upwind  profile.  Because  the  vertical  velocities,  w,  needed  to  determine  stream  lines  are 
incorrectly  assumed  to  be  zero  before  QUIC  performs  the  variational  procedure,  this  streamline 
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is  incorrectly  assumed  to  be  the  horizontal  dotted  line.  To  correct  this,  a  second  execution  of  the 
variational  procedure  is  performed.  After  the  second  variational  procedure,  the  streamline  is 
correctly  modeled  by  taking  into  account  the  non-zero  w  values.  Graphically,  this  is  shown  in 
Figure  1 1  with  the  modified  streamline  after  the  second  variational  procedure  shown  in  red.  It  is 
found  that  additional  runs  of  the  variational  procedure  produce  no  changes  in  the  origin  of  the 
tank  coordinate  system,  thus  the  variational  procedure  is  performed  only  twice.  The  procedure  is 
shown  in  Figure  12. 


Figure  11:  Illustration  of  the  results  from  the  second  (correctional)  variational  procedure  for  the  mixed 
layer  at  the  top  of  the  windbreak.  The  second  variational  procedure,  or  iteration,  accounts  for  vertical 
velocities  present  over  windbreaks,  and  as  a  consequence  more  accurately  places  the  origin  of  the  tank 
function  along  the  streamline  originating  at  the  windbreak  height  in  the  upwind  profile. 
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Location  of  streamline  located  at  Hnt,  far  upstream  and  the 
location  of  the  tanh  function’s  origin  on  that  streamline 


I 


Figure  12:  Flow  chart  of  procedure  utilized  to  simulate  a  windbreak  velocity  field, 
variational  procedure  is  performed  twice  in  QUIC. 


Note  that  the 


Turbulent  Field 


Similarity  solutions  for  windbreak  turbulence  have  been  explored  by  Judd  et  al.  (1996)  and 
Perera  (1981).  In  the  general  form  of  these  solutions,  the  mean  turbulence  quantities  are 
normalized  by  a  turbulent  velocity  scale,  Us.  The  vertical  dimension  is  normalized  by  23,  the 
thickness  of  the  shear  zone  downwind  of  the  top  of  the  windbreak,  and  the  origin  of  vertical  axis 
is  translated  from  the  ground  to  H.  Mathematically  this  is  expressed  as: 


where 
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(  =  (z~  H)/25. 


(9) 


Judd  et  al.  (1996)  took  t/|  to  be  the  magnitude  of  the  maximum  Reynolds  stresses,  it'w'  at  a 
given  downwind  distance  of  the  windbreak.  The  resulting  scaling  is  quite  good  and  shown  in 
Figure  13. 
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Figure  13:  The  turbulence  similarity  solution  proposed  by  Judd  et  al.  (1996). 

There  are  three  steps  to  parameterizing  the  turbulence  similarity  solution:  1)  producing  the 
mathematical  description  of  the  vertical  profile  given  in  Figure  13,2)  specifying  Us,  the  turbulent 
velocity  scale,  at  all  downwind  locations  and  3)  relating  all  turbulent  quantities,  it'2  and  w'2  to 
the  Reynolds  stresses,  it'w'. 

For  the  first  step,  the  functional  description  of  the  vertical  dependence  of  u'w' /U2  on  2  was 
chosen  to  be: 


fiOuw  =  Atanhtf  ~  c(A,B ))  +  BcoshT 2(^  -  c(A,B ))  +  C.  (10) 

Equation  10  includes  tank,  used  in  the  windbreak  mean  velocity  profile,  and  its  first  derivative, 
cosh'2.  The  constant  c(A,B)  ensures  the  maximum  of  Eq.  10  occurs  at  £,  =0.  Its  value,  c(A,B) 
=atanh(A/(2B)),  may  be  found  by  equating  the  vertical  derivative  of  Eq.  10  to  zero.  The 
parameters  A,  B,  and  C  subject  the  Eq.  10  profile  to  three  constraints. 

Atanh(oo)  +  B cosh-2  (oo)  +  C  =  u'w'upstream  (X  ->  oo)  (11) 

Atanh(—co)  +  Bcosh~2(—co)  +  C=  u'w' min  ->  —  oo)  (12) 

Atanh(c(A,  B))  +  Bcosh^^c^A,  B))  +  C  =  U2  (l;  =  0)  (13) 


Where  u'w' upstream  is  the  upstream  surface  layer  value  of  it'w'  and  u'w'min  is  the  Reynolds 
shear  stress  in  the  bleed  flow.  It  is  assumed:  u'w'min  —  apleed  flow  *  u'w'upstream.  Also  Eq.  10 
assumes  that  it'w'  is  a  small  non-zero  value  at  the  ground  instead  of  its  true  value,  zero.  The 
consequences  of  this  discrepancy  are  negligible  for  the  resolution  of  data  examined  in  this  work. 
Using  the  relations  tanh(c(A,B))  =  A/(2B),  1-tanh  =cosh'  and  c(A,B)  >  0  together  with  the  limits, 
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cosh~  (oo)  0,  cosh~  (-cc)  0,  tanh(co)  1 ,  tanh(-co)  -1  yield  the  values  for  the  constants 
which  are:  A  =  (1  -  u'w'min)/2  ;  C  =  (1  +  u'w'rnin)/2  ;  B  =  ((-C  +  Us )  +  ((C  -  Us  )2  - 
A2)05)/2. 


For  the  second  step,  parameterizing  Us,  the  turbulent  velocity  scale,  at  all  downwind  locations, 
x/H,  we  hypothesize  that  there  is  a  similarity  solution  U*s  —  f(x *).  A  model  of  this  solution  is 
shown  in  Figure  14.  U*s  and  x*  are  similarity  parameters  that  account  for  variation  of  windbreak 
aerodynamic  porosity  aobf  ,  and  windbreak  height  compared  to  upstream  roughness  ( H-d)/za .  l/s* 

is  the  difference  of  Us  and  its  undisturbed  upstream  value,  u'w'°UpStrearn,  nonnalized  by  the 
difference  of  Us  max  (the  maximum  value  of  Us  located  in  the  high  shear  zone  just  downwind  of 
the  windbreak  top)  and  w'w'ypStrearn.  We  can  then  write  this  nonnalized  scaling  velocity  as: 


JT *  Us—U  W  upstream 

Us  ~  j  ,0.5  5 

Us  max~u  w  upstream 

and  the  nonnalized  downwind  distance,  x*  as: 

x*  =  (Dao^(z0/(H  —  d))  75  +  £)  *  ( x/H  —  Fln((H  —  d)/z0)  +  G )  . 

The  variables:  D,  E,  F,  and  G  are  coefficients  to  be  detennined  by  a  least  squares  lit  to  data. 

This  parameterization  is  based  upon  the  data  sets  of:  Perera  (1981),  Bradley  and  Mulhearn 
(1983),  and  Judd  et  al.  (1996)  all  which  have  1)  of  the  relatively  high  spatial  resolution  of 
turbulence  measurements,  2)  the  span  in  range  of  windbreak  heights  nonnalized  by  upwind 
roughness  {H-d)/z0  =  30  to  ( H-d)/z0  =  600  3)  and  3)  the  large  range  of  optical  porosities  [i=0  to 
0=0.7.  Similar  relationships  have  been  developed  previously  (  Perera  1981)  but  are  valid  only  at 
distances  far  downstream  of  windbreaks,  x/H> 7.5  and  are  validated  with  a  single  dataset  instead 
of  multiple  data  sets  as  is  done  in  this  work.  The  result  of  the  least  squares  fit  of  x*and  is 
shown  in  Fig.  15. 
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Figure  14:  The  variables  x*  and  illustrating  the  recovery  of  Us,  the  turbulent  velocity  scale,  from  its 
maximum  value,  Us  max,  (occurring  on  the  downwind  side  of  the  windbreak  on  the  left)  to  its  upstream 

value,  u'w' uvstream  far  downwind  on  the  right.  The  recovery  is  divided  into  three  regimes:  the  first 
regime  exhibits  asymptotic  decay  from  U*s  =  1 ,  the  second  exhibits  rapid  decay  of  U*s ,  and  the  final 
exhibits  U*s  approaching  0  asymptotically. 
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1.2 


* 
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Figure  15:  The  calculation  of  the  turbulent  velocity  scale,  Us  _  utilizing  two  dimensionless  parameters: 
U s*  and  x*  .  x*  =  (ll.6a£jjy(z0/(//  -  d)>75  +  0.074)  *  (x/H  -  7.9 Zn((W  -  d)/z0)  +  21)  and  l/s*  = 

( Us/u'w'0UpStrearn  ~  1)/ (0.22  Zn((l  -  a0&/)2  *  z0/(H  -  d))  -  0.13),  where  aobf  is  the  aerodynamic 
porosity  of  the  bleed  flow  through  the  windbreak.  The  coefficients  are  obtained  from  least  square  fits  to 
data  from  Judd  et  al.  (1996),  Bradley  and  Mulheam  (1983),  and  Perrera  (1981).  t/Jand  x*are  related  by: 
U*s  =  -0.5  *  tanh(x*)  +0.5. 


Once  Us  has  been  parameterized,  the  mathematical  description  of  it'w'  is  complete  and  the 
remaining  elements  of  the  Reynolds  stress  tensor,  it'2  and  w'2,  can  also  be  parameterized.  The 
following  relations  are  used  by  Counihan  et  al.  (1974)  for  this  purpose  and  are  valid  for  the 
surface  layer: 


it'2  =  S.Ou'w', 
w'2  —  1.5u'w'. 

3.2.2.  Improved  Deposition  Modeling 

The  wind  tunnel  research  outlined  in  this  report  was  designed  to  help  develop  a  model  for 
particle  deposition  with  improved  physics  compared  to  current  models.  While  the  model  form  is 
presented  here,  the  basis  for  the  model  is  described  in  Section  3.1.1  in  the  results  section  of  the 
wind  tunnel  testing.  This  section  now  details  the  method  of  incorporating  the  resulting  model 
into  QUIC. 


(14) 

(15) 
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Particle  deposition  is  often  quantified  by  the  deposition  velocity;  Vd.  A  model  for  Vd  that  is 
commonly  utilized  is  that  of  Raupach  et  al.  (2001).  In  this  model,  a  particle  of  a  given  diameter, 
dp  and  traveling  at  a  velocity  of  U  interacts  with  a  vegetative  element  of  size  de.  Vegetative 
orientation  is  neglected.  The  Stokes  number,  Stk,  for  the  particle  element  pair  is  calculated  as: 

Stk  =  (ppdp)/(18pava)[de/2I7],  (16) 

•5 

where  pp  is  the  density  of  the  particle  (assumed  1000  kg/m  ),  pa  is  the  density  of  the  air,  and  va  is 
the  kinematic  viscosity  of  the  air.  Once  St  is  calculated,  Vd  can  be  calculated  by 

DF=(Stk/ (St  +  p))q  (17) 

Vd  —  U  ■  DF,  (18) 

where  DF  is  the  deposition  fraction  originally  called  deposition  efficiency  by  Raupach  et  al. 
(2001),  p  and  q  are  constants  and  Raupach  et  al.  (2001)  utilized  0.8  and  2  respectively.  In  a 
Lagrangian  reference  frame  the  concentration,  C,  decreases  as: 

dc/dt  —  —  VdyC.  (19) 

where  y  is  the  vegetative  surface  density  (vegetative  surface  area/unit  volume  of  vegetation)  his 
can  be  integrated  along  the  particle  path  assuming  Vd, ,  and  U  to  be  constant.  This  yields: 

er  =  Cx/C0  =  exp((—  VdySb/U),  (20) 

Where  5),  is  the  distance  traveled  by  the  particle  within  the  vegetative  canopy  during  the  time  of 
integration  and  Cj  and  C0  are  the  concentrations  at  the  beginning  and  end  of  the  period  of 
integration  respectively.  The  details  of  this  model  are  presented  in  Amatul  (2006). 

An  important  deliverable  from  the  wind-tunnel  research  was  an  improved  parameterization  for 
deposition  fraction,  DF,  that  accounts  for  deposition  enhancement  by  turbulence  (see  results  in 
section  3.1).  DF  is  a  function  of  the  turbulence  and  the  Taylor  micro-scale  (k),  which  is  an 
important  parameter  for  quantifying  the  turbulence  occurring  within  vegetation.  It  is  very 
difficult  to  directly  compute  a  general  Taylor  micro-scale.  However,  if  one  assumes  that  the 
turbulence  is  locally  isotropic,  following  Tennekes  and  Lumley  (1972)  k  can  be  calculated  within 
QUIC’s  existing  turbulence  framework  via  the  following  relationship: 


1 


In  Eq.  21,  s  is  the  dissipation  rate  of  turbulent  kinetic  energy  and  uj  is  the  squared  local 
turbulence  intensity. 

The  final  step  is  to  use  the  wind-tunnel  developed  expression  (Eq.  47)  for  DF  instead  of  the 
traditional  Eq.  17.  This  is  described  below. 
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3.3.  Full-scale  Deposition  Field  Experiments 

To  better  understand  transport  and  deposition  in  vegetative  canopies  and  windbreak  two  field 
experiments  were  conducted.  The  first  was  at  Hanford,  Washington  (WA).  This  site  featured 
native  vegetation  in  fonn  of  a  canopy.  The  second  utilized  the  Raft  River  Windbreak  near  Malta, 
Idaho  (ID).  Both  experiments  included  high-resolution  measurements  of  the  flow  field  as  well  as 
concentrations  fields  of  particulate  matter. 

3.3.1.  Hanford,  Washington  Semi-arid  Vegetation  Experiment 

A  field  study  was  performed  at  the  Pacific  Northwest  National  Laboratories  Hanford  site 
(46°34’3.53"N  and  1 19°36’5.44"W,  -220  m  above  sea  level).  The  dates  of  the  campaign  were 
from  June  6  to  June  13,  2011.  During  the  field  campaign,  two  Intensive  Observation  Periods 
(IOPs)  were  conducted  with  simulated  road  traffic.  Table  2  summarizes  the  IOPs.  The  study  was 
conducted  in  conjunction  with  J.  Gillies  and  V.  Etyemezian’s  group  from  the  Desert  Research 
Institute  (RC-1729).  The  region  around  the  measurement  site  contained  native  vegetation  in  all 
directions  of  about  1 .4  m  in  height.  The  road  was  oriented  due  east- west  and  was  unpaved. 

The  two  objectives  for  the  field  study  were: 

1)  Objective  1:  To  verify  the  ability  of  QUIC  to  estimate  PMiotransport  and  deposition 
within  a  vegetative  canopy  typical  of  many  military  facilities.  The  QUIC  simulation 
was  perfonned  in  real  time  during  an  IOP  on  June  7,  2013  where  the  aerosol  sensors 
were  located  based  on  the  results  of  the  QUIC  simulation  of  the  PMio  dispersion  and 
removal  at  the  site  in  real-time.  This  represented  a  demonstration  of  the  prototype 
field  model  for  planning  and  prediction. 

2)  Objective  2:  To  better  understand  the  behavior  with  PMio  within  vegetative  canopies. 

Fugitive  dust  was  generated  by  vehicular  traffic,  a  large  van  or  truck,  making  discrete  separate 
trips  on  the  unpaved  road.  In  between  trips,  time  was  allowed  to  pass  so  the  dust  cloud  resulting 
from  one  trip  could  clear  all  the  equipment  before  the  cloud  resulting  from  the  next  trip  entered 
the  equipment  array. 

Supporting  Activities 

To  relate  the  fugitive  emissions  to  models  and  other  field  studies,  a  number  of  site  dependent 
measurements  were  necessary.  These  included  site-specific  vegetative  measurements  and 
meteorology  measurements.  The  methodology  utilized  to  document  the  roughness  and 
meteorology  conditions,  are  now  presented. 


IOP 

Purpose 

IOP  Date/  Time 

Atmospheric  Conditions 

zref=  3.9  m 
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1 

Use  QUIC  to  design  field 
study  and  make  real-time 
sensor  position  adaption 

June  7,  201 1/12:50  to  14:50 
Pacific  Daylight  Savings  time 
(PDT). 

uref  =  5.4  m  s'1 

L=  -47.4  m 
u*  =  0.61  m  s'1 
tke  =  0.4  nf  s'2 
z0  =  0.07  m 

2 

Validate  QUIC  with  high 
resolution  canopy  PMi0 
profile 

June  13,2011/  16:30  to  17:30 
PDT. 

uref=  6.8  m  s'1 

T=  - 1 1 7  m 
u*  =  0.71  m  s'1 
tke  =  0.3  nr  s'2 
z0  =  0.07  m 

Table  2:  Summary  10P  table  for  the  Handford,  WA  experiment.  Averages  were  taken  over  the  10P 
duration. 

Meteorology 

For  the  field  study,  it  was  important  to  carefully  characterize  the  basic  transport  mechanisms,  in 
particular:  horizontal  advection  and  vertical  transport.  To  characterize  these  processes 

meteorological  data  were  taken  within  and  above  the  vegetative  canopy  at  one  tower  using  four 
Campbell  Scientific  3D  anemometers  (CSTA3).  The  3D  sonic  anemometers  measured  the  three 
velocity  components  and  virtual  temperature  at  a  rate  of  10  Hz.  A  Campbell  Scientific  CR5000 
data  logger  was  utilized  to  record  and  process  the  data.  In  order  to  ensure  reliability  of  the 
measurements,  all  of  the  anemometers  and  data  loggers  were  within  the  recommended 
manufacture’s  calibration  period.  The  tower  location  and  anemometer  heights  are  noted  in  Table 
3.  The  anemometers  are  shown  in  Figure  16.  The  location  of  the  sonic  anemometers,  37  m 
downwind  of  the  road,  was  chosen  to  allow  the  airflow  to  recover  from  the  transition  of  the 
exposed  ground,  characteristic  of  the  road,  to  the  vegetation,  which  was  about  1.4  m  in  height. 
The  data  utilized  in  making  this  choice  are  found  in  Arya  (2001). 

The  sonic  anemometers  make  possible  a  characterization  of  the  transport  mechanisms  within  and 
above  the  vegetation.  Parameters  needed  to  perform  QUIC  simulations  were  calculated  from  the 
sonic  anemometer  data  in  real  time  during  the  experiment.  This  enabled  fast  computer 
simulations  to  optimize  sensor  placement  for  IOPs  involving  Objective  1. 


Equipment 
(Number  of  Sensors) 


Downwind  Distance  (m) 


Sensor  Heights  (m) 
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3D  Sonic  Anemometers  (4) 


37.0 


0.68,  1.40,2.00,3.90 


Table  3:  Location  of  anemometers  that  characterized  meteorological  conditions  during  the  Handford  WA 
field  study. 


Figure  16:  Photo  of  the  3D  Sonic  Anemometer  Tower  at  the  Hanford  field  study  site. 


Site  Vegetation 
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The  density  of  vegetation  was  characterized  using  the  Leaf  Area  Index  (LAI).  LAI  is  the  ratio  of 
plant  surface  area  to  ground  surface  area  and  is  dimensionless  (Breda,  2003).  The  current  project 
utilized  two  methods  for  estimating  LAI  at  the  field  site.  The  first  was  remote  sensing  based, 
utilizing  the  United  States  Geological  Survey  (USGS)  MODIS  database.  This  is  located  online  at 
the  URL  of: 

http://daac.ornl.gov/cgibin/MODIS/GLBVIZ  1  Glb/modis  subset  order  global  col5.pl.  (Last 
accessed  in  2011)  The  second  method  was  an  indirect  optical  method  to  detennine  LAI  using  a 
Licor  LAI-2000®. 

In-field  calibration  of  the  LAI-2000  included  a  diffuse  shortwave  radiation  measurement  to  be 
made  at  sunrise,  sunset,  or  on  a  cloudy  day.  The  vertical  profile  of  LAI  was  measured  at 
multiple  locations.  This  was  performed  by  dividing  the  height  of  the  vegetation  canopy  into  four 
approximately  equal  vertical  intervals  defined  by  the  points:  0 Hc,  0.25 Hc,  0.5 Hc,  0.75//c.  An 
advantage  of  these  onsite  measurements  was  that  the  distribution  of  LAI,  or  vegetative  density 
with  height  could  be  detennined.  Using  the  MODIS  satellite  imagery  database  method,  the 
distribution  with  height  was  unknown.  The  variation  of  LAI  with  height  may  be  important  for 
simulation  purposes.  These  LAI  profile  measurements  were  performed  at  locations  along  the  x- 
axis  in  Figure  19  in  2  m  intervals  from  the  roadside  to  the  distance  of  the  farthest  downwind 
measurement  (~60  m  downwind  of  the  road).  Also,  a  number  of  measurements  were  performed 
parallel  to  the  road  direction  at  37  m  downwind  of  the  road  the  location  to  ensure  that  the  LAI 
calculations  include  vegetation  that  may  have  influenced  fugitive  dust  emissions  if  the  wind  was 
not  perfectly  perpendicular  to  the  road  direction.  In  total,  48  measurements  were  made.  The  total 
distance  from  the  sonic  tower  the  measurements  covered  was  dictated  by  the  experimental  results 
of  Cheng-I  et  al.  (2000).  This  distance  included  all  the  vegetation  that  influenced  the  sonic 
measurements. 


Fugitive  Dust  Concentration  Measurements 

Measurements  of  dust  concentrations  were  made  using  TSI  DustTrak®  8520  Aerosol  Particulate 
Monitor  optical  sensing  devices  sampling  at  1  Hz.  The  DustTrak  is  a  rugged  portable  device  that 
uses  light  scattering  of  dust  to  measure  concentrations.  It  is  capable  of  measuring  concentrations 
ranging  from  0.001  to  100  mg/m3  and  has  been  utilized  in  previous  studies  of  road  dust 
emissions  (Veranth  et  al.,  2003).  PMio,  inlets  were  utilized  in  this  experiment.  The  DustTrak 
sensors  were  rented  from  Galson  Labs  (http://galsonlabs.com/)  and  Argus-Hazco  (http://argus- 
hazco.com/).  The  companies  calibrate  the  instruments  prior  to  shipment.  The  standard 
calibration  procedure  used  to  calibrate  the  instruments  is  found  in  the  TSI  8520  DustTrak 
Manual.  Further,  colocation  tests  were  performed.  The  data  from  are  shown  in  Appendix  A. 
Typical  inter-instrument  variability  characterized  by  the  standard  deviation  of  mean 
concentration  measurements  was  0.06  mg/m  compared  to  typical  maximum  concentrations  of 
40  mg/m3. 

Objective  1 

As  previously  noted,  QUIC  was  previously  capable  of  modeling  the  dispersion  of  PMio  that 
occurs  in  vegetative  canopies  and  has  a  model  for  vegetative  deposition  (Amatul  2006).  The 
dispersion  parameterizations  are  an  enhancement  of  the  work  of  Cionco  (1965)  and  utilize  his 
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vertical  velocity  profile  for  the  removal  of  momentum  from  the  wind  field.  This  new  vertical 
profile  is  used  in  conjunction  with  a  mixing  length  turbulence  model  to  estimate  dispersion. 

In  field  studies,  the  dispersion  of  PMio  was  measured  utilizing  a  flux  methodology.  Horizontal 
fluxes  of  the  PMio  component  of  fugitive  dust,  denoted  by  F,  are  common  in  estimating  the 
transport  of  fugitive  dust  (Veranth  et.  al  2003;  Etyemezian  et  al.  2004;  Speckart  et  al.  2013). 
Differences  in  F  at  two  points  are  equal  to  the  removal  of  the  PMio  between  the  two  points.  F  is 
the  integrated  product  of  the  local  concentration  and  horizontal  wind  speed.  The  integration  is 
perfonned  from  the  ground  to  the  cloud  top,  and  then  summed  over  the  time  of  travel  for  the 
plume  to  travel  past  the  downwind  PMio  sensors.  Both  wind  speed  and  concentration  have  been 
measured  experimentally  during  all  the  IOPs  for  this  current  work.  Owing  to  the  discrete  nature 
of  the  measurements  in  section  field  studies,  modeling  assumptions  are  needed  to  estimate  fluxes. 
This  project  utilizes  the  method  proposed  by  Speckart  et  al.  (2013)  (This  method  is  described 
later  in  this  section.)  All  of  the  methods  yield  the  same  result  when  a  sufficiently  large  number  of 
sensors  are  utilized.  When  fewer  sensors  are  used,  the  proposed  method  is  more  robust  than 
previously  utilized  methods. 

In  the  absence  of  sources,  the  horizontal  flux  will  be  a  non-increasing  function  with  increasing 
downwind  distance  from  a  road.  Mechanisms  such  as  vegetative  deposition  and  gravitational 
settling  can  remove  fugitive  dust  and  consequently,  reduce  horizontal  fluxes.  This  reduction  can 
be  modeled  in  QUIC  as  illustrated  in  Figure  17. 


Downwind  Distance  (m)  (Along  x  axis  in  Figure  18) 

Figure  17:  Illustration  of  a  QUIC  Dispersion  model  simulation  of  the  reduction  of  roadside  horizontal 
flux  F/Froadside  for  the  Hanford  experiment. 

To  facilitate  field  and  simulation  comparison,  the  horizontal  flux  is  generally  normalized  by  the 
roadside  flux,  Fma(islde-  As  one  may  deduce  from  Figure  17,  the  majority  of  the  removal  occurs 
near  the  road  where  the  majority  of  the  dust  cloud  interacts  with  vegetation  and  the  ground.  As 
downwind  distance  increases,  the  removal  mechanisms  become  less  important  and  the  flux 
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becomes  nearly  constant  at  a  value  denoted  as  Fmin.  Owing  to  this  behavior,  a  new  quantity,  the 
fraction  removal,  Frac,  can  be  defined  as: 

Frac(x)  =  F^~F^  .  (22) 

^roadside  ^  mm 


Frac  has  a  maximum  of  1  at  the  roadside  and  a  minimum  of  0  far  downwind  from  the  road.  The 
behavior  of  this  fraction  of  removal  is  illustrated  in  Figure  18.  Frac  may  be  physically 
interpreted  as  the  remaining  fraction  of  removal  that  will  occur  as  one  travels  downwind.  At  the 
roadside,  all  removal  mechanisms  are  located  downwind,  hence  Frac- 1.  At  a  position  where 
Frac= 0.5,  half  of  the  removal  is  located  upwind  and  half  is  located  downwind,  and  when  Frac=0, 
all  of  the  removal  mechanisms  are  located  upwind. 


Downwind  Distance  (m)  (Along  x  axis  of  figure  5) 


Figure  18:  Reduction  in  normalized  dust  concentration  (Frac)  with  increasing  downwind  distance  from  a 
QU1C  simulation.  DustTrak  tower  locations  for  the  Hanford  experiment  were  determined  in  real-time 
using  the  Frac  curve  as  predicted  by  the  QUIC.  Values  of  1,  0.5,  0.1,  were  chosen  for  comparison 
between  QUIC  simulation  and  the  field  study.  These  values  indicate  the  fraction  of  the  total  removal  that 
is  predicted  to  occur  as  the  cloud  advects  further  downwind.  The  total  removal  may  be  significant  or 
insignificant,  Frac,  indicates  the  distribution  of  the  total  removal  with  downwind  distance. 


The  placement  of  the  three  DustTrak  towers  during  the  experiment  was  determined  using  QUIC 
simulations  such  that  the  first  tower  was  placed  roadside  where  Frac  =  1,  the  second  tower 
where  Frac= 0.50,  and  the  third  tower  was  located  where  Frac  =  0.10.  The  experimental  set-up 
is  summarized  in  Figure  19.  The  physical  locations  for  the  placement  of  the  sensors  during  the 
IOPs  is  shown  in  Table  4  along  with  the  date,  number  of  runs,  and  atmospheric  stability. 
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V 


x  axis 


Figure  19:  Schematic  of  equipment  deployment  for  Objective  1  lOPs.  The  DustTrak  towers  are 
dimensioned  by  the  placement  criteria,  DustTrak  towers  located  at  the  physical  locations  where  F rac  =1, 
0.5,  and  0.1  according  to  QU1C  simulations. 


IOP  DATE 

Atmospheric 

Stability 

Number  of  Vehicle 
Runs 

Downwind  Distance 
(m)  (Along  x  axis  of 
Fig.  1)  of  DustTrak 
monitor  towers.  All 
towers  employ 
sensors  at  heights  of 
0.58,  1.20,1.70  m. 

June  7,  201 1 

Unstable 

32 

10.0,  16.4,  66.4 

Table  4:  The  placement  of  DustTrak  Sensors  : 

or  objective  1  IOP  during  the  Hanford  canopy  c 

study.  Note  that  the  downwind  DustTrak  tower  location,  right  most  column,  moves  according  to  QU1C 
simulation  results. 


QUIC  dispersion  and  deposition  simulations  were  run  in  order  to  field  test  the  QUIC  model  and 
optimize  field  sensor  locations  using  real  time  monitored  variables  from  sonic  anemometer.  The 
needed  data  included:  horizontal  mean  wind  at  the  four  heights,  the  friction  velocity  (u*),  the 

absolute  near  surface  temperature  Ta,  the  kinematic  heat  fluxes  Qo  —  w' (computed  from  the 
sonic  temperature  and  vertical  velocity  correlation),  the  Monin-Obukhov  length  scale 
L  =  ull[KQ0(glTc)\  vegetative  height,  and  LAI.  These  meteorological  variables  are  described 
in  Arya  (2001).  With  the  experiment  configuration  utilized  in  Hanford,  these  variables  were 
obtainable  just  before  the  first  vehicle  trip  over  the  unpaved  road.  The  meteorological  data  were 
based  on  5  min  averages. 
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To  estimate  the  robustness  of  the  QUIC  simulation  results,  Eq.  22  can  be  solved  for  Fmin.  This 
yields: 


F 

min 


FracF  —  F 

1  rc lL'1  roadside  1 

Frac  —  l 


(23) 


An  estimate  for  Fmin  was  computed  using  the  IOP  data  obtained  from  the  second  and  third  towers. 
Another  estimate  for  Fmin  was  obtained  from  the  QUIC  simulation  of  the  IOP,  with  identical 
meteorological  and  roughness  conditions.  The  estimated  accuracy  of  the  QUIC  model  increases 
as  the  difference  of  the  two  estimates  for  Fmin  decrease. 

Objective  2 

Past  studies  have  focused  on  fugitive  dust  emissions  by  examining  PMio  concentrations  at  all 
heights,  from  the  ground  to  the  cloud  top  (Veranth  et  al.  2003;  Etyemezian  et  al.  2004;  Speckart 
and  Pardyjak  2011).  As  a  consequence  of  limited  resources,  a  limited  number  of  sensors  have 
been  used  to  measure  vertical  profiles  of  PMio  concentrations  within  vegetative  canopies. 
Models  for  fugitive  dust  concentrations  resulting  from  emissions  have  been  implemented  above 
vegetative  canopies  with  success,  such  as  exponential  and  step  profiles  (Veranth  et  al.  2003)  but 
may  be  inappropriate  within  the  canopies.  Issues  such  as  the  boundary  condition  for  PMio 
concentration  at  the  ground  are  important  for  simulation  but  to-date  are  relatively  poorly 
understood. 

Concentrations  within  vegetative  canopies  are  generally  assumed  to  decrease  with  increasing 
height  for  near  ground  sources  (Veranth  et  al.  2003).  This  is  a  common  assumption  in  the 
calculation  of  F  (e.g.  Veranth  et  al.  (2003)).  This  is  the  case  despite  that  Lagrangian  dispersion 
models  that  suggest  otherwise  (Aylor  and  Flesch  2001).  Currently,  there  is  a  lack  of  high- 
resolution  concentration  data  measured  in  the  field  (e.g.  more  than  3  points  to  estimate 
concentration  profiles  within  canopies).  Hence,  we  designed  an  experiment  to  measure  the 
vertical  profile  of  PMio  with  a  tower  consisting  of  nine  DustTrak  monitors.  Their  heights  and 
downwind  distance  are  shown  in  Table  5. 


Date 

Number 
of  Vehicle 
Runs 

Downwind 
Distance  (m) 

Sensor  Heights  (m) 

June  13, 2011 

50 

32 

0.04,0.28,0.58,0.71,0.89,  1.17, 
1.48,  1.78 

Table  5:  3D  Sonic  and  DustTrak  Sensor  placement  for  the  canopy  concentration,  objective  2  IOP 
measured  during  16:00  to  18:00  PDT  and  featuring  unstable  atmospheric  conditions. 

A  photograph  of  the  tower  with  DustTraks  is  shown  in  Figure  20.  Utilizing  the  DustTrak 
concentration  data,  a  detailed  vertical  concentration  profile  can  be  produced. 
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Figure  20:  Photograph  of  DustTraks  during  the  high-resolution  vertical  profile  concentration,  objective 
two  IOP  on  June  13,  2012  at  Handford,  WA. 


Method  for  Calculating  Horizontal  Flux,  F  for  Handford,  WA  Field  Study 

PMio  concentration  and  wind  data  can  be  used  in  multiple  ways  to  calculate  the  horizontal  flux. 
Horizontal  flux  at  x  was  calculated  as: 

Fx  —  fZ_Q  f^max  C(x,  z,  t)it(x,  z,  t)dtdz.  (24) 

tmax  is  the  time  for  the  plume  to  travel  more  than  100  m  down  wind.  In  the  subsequent  discussion, 
cint(x,  z)  —  ft=Qmax  C(x,z,t)dt,  and  is  utilized  to  apply  a  continuous  dispersion,  such  as  a 
Gaussian  Plume  model,  to  data  obtained  from  discrete  vehicle  trips.  The  new  method 
incorporates  additional  physics  and  complexity  to  substitute  for  large  amounts  of  concentration 
data.  In  addition  to  wind  and  PMio  concentration  data,  it  requires  an  estimate  of  atmospheric 
stability  through  the  Monin-Obukhov  length  scale,  L  (Ayra  2001).  L  can  be  calculated  using 
eddy  covariance  methods  (see  section  2.3.1),  profiles  of  mean  velocity  and  temperature  (Ayra, 
2001)  or  the  Pasquill  stability  chart  (Seinfield  and  Pandis,  1998).  For  the  work  presented  here, 
the  eddy  covariance  method  was  used.  The  new  method  utilized  the  Gaussian-plume  model 
proposed  by  Van  Olden  (1978): 

cint(x,  z)/Q  =  (A(x)  exp  [(—5  z/z  (x))s]  ),  (25) 

where  Q  is  the  source  strength  (kg  m'1),  A(x)  is  a  coefficient  quantifying  the  downwind  decrease 
of  cint(x,  z),  s,  the  shape  factor,  was  a  constant  that  varied  from  1  to  2  (Venkatram  2004).  If  5=1 
the  profile  is  exponential,  if  5=2  the  profile  is  Gaussian.  Venkatram  (2004)  assumed  5=1.5  as 
suggested  by  Van  Olden  and  this  work  utilized  that  assumption.  B  is  a  constant  defined  as: 

B  =  r(2/s)/r(l/s),  (26) 

where  f  was  the  gamma  function. 
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The  height,  z(x),  is  the  vertical  center  of  mass  of  the  PMio  plume.  Its  definition  is: 


z(x)  = 


/0°°zcmt(x,  z)dz 

/q  C int  (x,  z)  dz 


(27) 


Eq.  27  suggests  that  the  vertical  profile  of  ln(cint(x,  z))  and  zs  is  linear  and  consequently  its 
characterization  is  simpler  than  the  non-linear  profile  of  cint(x,z)  and  z.  As  a  result,  the  new 
method  estimates  the  profile  of  Cjnt(x,  z)  and  z  indirectly  by  first  calculating  the  profile  of 
ln(cjnt(x,  z))  and  zs  and  subsequently  computing  the  profile  of  cint(x,z)  and  z  by  taking 
exponentials  and  roots  of  the  respective  variables. 


The  new  method  divides  the  integrated  concentration  profile  into  three  zones:  1)  the  zone  below 
the  lowest  DustTrak  2)  the  zone  between  the  lowest  and  highest  DustTrak  monitor  3)  the  zone 
above  the  highest  DustTrak  monitor.  This  is  shown  in  Figure  21.  For  the  lowest  zone, 
In  (cjnt(x,  z))  is  assumed  to  be  constant  from  the  ground  to  the  bottom  DustTrak  monitor.  For 
the  intennediate  zone,  within  the  DustTrak  monitors,  ln(cjnt(x,  z))  is  assumed  to  vary  linearly 
as  a  function  of  z?  between  adjacent  DustTrak  monitors.  Above  the  topmost  DustTrak  monitor, 
ln(cint)  is  assumed  to  decrease  linearly  with  increasing  z\ 


Linear  decrease  of  In  cint  with  increasing  zs  above  top 
DustTrak  monitor  as  indicated  by  the  Van  Olden  (1978) 
model. 


Figure  21:  Conceptual  illustration  of  the  methodology  used  to  calculate  the  horizontal  flux,  F  from 
DustTrak  data  .  £ 


Above  the  top  DustTrak  monitor,  we  had  to  calculate  the  rate  of  decrease  of  ln(c,nt)  with 
increasing  zs  (denoted  as:  din (cint)/dzs).  For  a  given  downwind  distance,  din (cint)/dzs  is  a 
constant.  In  other  words,  the  vertical  decrease  of  ln(Cjnt)  with  increasing  zs  is  linear. 
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d  ln(cint) / d zs  can  be  calculated  by  taking  logarithms  and  then  derivatives  of  both  sides  of  Eq. 
25. 


d  In  cint/dzs  =  —  ■ 


(28) 


After  deriving  Eq.  28  it  remained  to  calculate  z(x)  .  The  downwind  dependence  of  z(x)  is 
modeled  by  the  ordinary  differential  equation  (ODE)  (Venkatram  2004): 


dz 

dx 


K(qz)/(u(qz)qz), 


(29) 


where  K  was  the  vertical  turbulent  diffusion  coefficient  and  q  is  defined  by: 


q  =  [s{r(2/s)/r(i/s)}s]*/n-*). 


(30) 


In  Eq.  28,  q  was  a  constant  that  relates  the  vertical  dispersion  to  distance  traveled  downwind. 
The  initial  condition  assumed  for  Eq.  29  was  z(0)  =  lm.  Vertical  profiles  for  K(qz)  and 
fi(qz)  are  needed  as  inputs  for  Eq.  29.  Venkatram  (2004)  inserted,  K (qz)  =  KU,qz,  for  a  neutral 
surface  layer  and  a  power  law  wind  profile  for  fi(z)  yielding: 

r  u  „i1/(p+ i) 

^  =  (1  /q)  [(p  +  l)xq  Tr-xZrJ  .  (31) 

This  was  extended  to  stable  atmospheric  conditions  using,  K(qz)  —  (/cu*qz)/(l  +  4.7 (q  z/L)): 


zp+ 1  4.7  qz^+2  KU^ 
p  +  1  +  (p  +  2)L  urqP 


(32) 


Eq.  32  can  be  solved  by  a  simple  bisection  root  finding  technique  (Chapra  and  Canale,  2006). 
For  unstable  conditions,  K(qz )  =  (/ai*qz)/(l  —  15 (qz/L))1^4,  would  have  yielded: 


dz 

dx 


p 

KU^Zy 

urqP 


(l  -  15(qz/L))1/4. 


(33) 


Eq.  32  can  be  solved  with  a  typical  integration  technique  such  as  a  fourth  order  Runge-Kutta 
method  (Chapra  and  Canale,  2006). 

The  vertical  profiles  of  cint(z )  were  obtained  by  taking  the  exponent  of  Zn(cjnt(z))  and  the  root 
of  zs . 
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3.3.2.  Raft  River  Windbreak  Experiment  -  Malta,  Idaho 

To  validate  the  new  models  implemented  in  QUIC  a  field  study  was  perfonned  during  July  of 
2012  at  the  Brigham  Young  University  Skaggs  Research  Ranch,  Cassia  County,  Idaho  located  in 
the  Raft  River  valley  of  southern  Idaho  (42°30’05.28"  N  1 13°24’29.92"  W)  near  Malta  ID.  This 
site  was  selected  because  in  1994  the  United  States  Department  of  Agriculture  (USDA)  Soil 
Conservation  Service,  Brigham  Young  University,  and  the  Idaho  Department  of  Transportation 
partnered  together  to  plant  several  sections  of  living  snow  fences  near  the  intersection  of  Idaho 
Hwy  81  and  US  Interstate  84  to  mitigate  the  effects  of  dust  storms  in  the  summer  and  snow  drifts 
in  the  winter  on  the  nearby  highways.  The  windbreaks  are  composed  of  two  sections  of  two  rows 
each  of  Rocky  Mountain  Juniper  ( Juniperus  scopulorum )  and  Siberian  Pea  ( Caragana 
arborescens ).  Supported  by  irrigation,  the  living  snow  fences  have  grown  to  -4.0  m  in  height  as 
shown  in  Figure  23. 

The  same  types  of  instruments  were  from  measuring  wind  speed  and  PMio  as  described  above 
for  Handford.  Again,  all  CSAT3  anemometers  were  within  the  manufacturers  recommend 
calibration  period.  The  DustTrak  sensors  were  rented  from  Galson  Labs  (http://galsonlabs.com/) 
and  Argus-Hazco  (http://argus-hazco.com/).  The  companies  calibrate  the  instruments  prior  to 
shipment.  The  standard  calibration  procedure  used  to  calibrate  the  instruments  is  found  in  the 
TSI  8520  DustTrak  Manual.  Further,  colocation  tests  were  performed.  The  data  from  are  shown 
in  Appendix  A.  Typical  inter-instrument  variability  characterized  by  the  standard  deviation  of 
mean  concentration  measurements  was  0.06  mg/m3  compared  to  typical  maximum 
concentrations  of  40  mg/m  .  Wind  data  were  taken  from  21:15  to  22:30  MDT  during  the  evening 
of  July  12,  2012.  Due  to  equipment  malfunction,  PMio  data  was  not  collected  during  this  period 
and  only  turbulence  data  were  recorded.  PMio  data  were  taken  on  July  27,  2012  from  12:00  to 
17:00  MDT. 

To  measure  the  mean  wind  and  turbulence  field  upstream,  within,  and  downstream  of  the 
windbreak,  five  towers  were  erected  with  Campbell  Scientific  CSAT3  sonic  anemometers 
(Campbell  Scientific,  Logan,  Utah).  Tower  A  was  placed  25  m  upstream  to  measure  free-stream 
turbulence  and  Tower  C  was  placed  in  between  the  two  windbreaks  to  measure  this  sheltered 
area.  Towers  D,  E,  and  F  were  placed  successively  downwind  of  the  windbreak  to  measure  the 
start  of  flow  recovery  to  upstream  conditions.  Table  6  gives  specific  coordinates  for  each  sonic 
anemometer,  DustTrak  and  windbreak.  Tower  B  had  no  3D  sonic  anemometers. 

PMio  data  were  taken  on  all  towers  except  tower  F.  Locations  are  documented  in  Figure  22  and 
Table  6.  The  objective  in  choosing  the  sensor  locations  was  to  characterize  as  much  as  possible 
the  approach  flow,  shelter  zone,  displaced  flow,  and  recovery  zones  that  are  noted  in  Figure  8. 
Both  fields  provide  verification  to  the  flow  models  developed  and  implemented  in  QUIC  for  this 
project. 
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Figure  22:  The  configuration  of  3D  sonic  anemometers  (Sonic)  and  DustTraks  for  the  Malta,  ID  Raft 
River  Windbreak  field  study. 
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Tower  A  (-34.1  m) 

Heights 

Sonic 

3.25  m 

DT 

3.15  m 

Tower  B  (1  m) 

Heights 

DT 

0.78  m 

DT 

1.39  m 

Tower  C  (13.8  m) 

Heights 

DT 

0.8  m 

Sonic 

1.03  m 

DT 

3.21  m 

Sonic 

3.31  m 

Tower  D  (27.8  m) 

Heights 

DT 

1.1  m 

Sonic 

1.25  m 

DT 

3.1  m 

Sonic 

3.25  m 

DT 

5.5  m 

Sonic 

5.6  m 

DT 

8.55  m 

Sonic 

8.7  m 

Tower  E  (41.1  m) 

Heights 

DT 

1.17m 

Sonic 

1.09  m 

DT 

3.08  m 

Sonic 

3.28  m 

Tower  F  (55.4  m) 

Heights 

Sonic 

3.25  m 

Table  6:  Location  of  equipment  deployed  for  the  Raft  River  Windbreak  field  study.  Experiment 
schematic  shown  in  Figure  22. 


Figure  23:  Photos  of  the  Raft  River  Windbreak  at  the  field  site  near  Malta  ID:  (left)  upstream  side  with 
unpaved  road  shown  and  (right)  downstream  side  of  the  windbreak. 
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3.4.  Development  of  Simple  Empirical  Formula  for  Native  Vegetative  Removal 

Prior  to  the  start  of  this  project,  Pace  (2005)  developed  the  simple  conceptual  model  of  the 
removal  of  dust  by  vegetation  shown  in  Figure  24.  This  qualitative  model  indicates  that  increases 
in  native  vegetative  height  and  thickness  increase  PMio  removal.  In  this  model  captured  fraction 
and  transmitted  fractions,  ( CF  and  TF  respectively  TF  =  (1  -CF))  of  vehicle  generated  PMio 
traveling  from  the  road  to  100  m  downwind  depends  upon  native  vegetative  height  and  thickness. 
CF  increases  with  increasing  height  and  density  of  downwind  vegetation.  In  this  section,  we 
present  a  simple  model  using  our  simulation  and  field  results  that  provides  quantitative  estimates 
consistent  with  Fig.  24,  and  can  be  easily  used  for  planning  purposes. 


Figure  24:  Recapitulation  of  the  conceptual  model  of  Pace  (2005)  that  indicates  that  the  transmitted 
traction  of  a  vehicle-generated  plume  is  a  function  of  the  height  and  density  of  downwind  roughness 
elements. 

We  utilized  a  simple  canopy  model  for  the  native  vegetation  in  which  the  canopy  height  is 
horizontally  homogeneous,  and  vegetative  density  is  assumed  constant  with  height.  For  such  a 
canopy  TF  and  CF,  can  be  calculated  by  the  following  simple  equation: 

TF  =  (1  —  CF)  =  IF  *  (1  —  canopy  removal  efficiency)  +  (1  —  IF)  (34) 
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Where,  IF  is  the  fraction  of  the  initial  dust  cloud  that  interacts  with  the  canopy  (the  interactive 
fraction)  and  the  canopy  removal  efficiency  defined  as  the  fraction  of  IF  that  is  removed  by  the 
canopy.  The  first  term  on  the  far  right  hand  side  (RHS)  is  the  fraction  of  the  initial  plume  that 
interacts  with  the  canopy  and  is  subject  to  removal  and  the  second  term  on  the  RHS  indicates  the 
fraction  of  the  initial  plume  that  is  not  influenced  by  the  canopy.  As  the  canopy  height  is 
increased  a  larger  and  larger  fraction  of  the  roadside  plume  interacts  with  the  canopy,  i.e.  the  first 
tenn  on  the  RHS  increases.  As  the  canopy  becomes  denser,  the  canopy  removal  efficiency 
generally  becomes  higher. 

To  quantify  Eq.  34,  two  canopy  parameters,  H*  and  T*,  were  introduced  in  Pardyjak  et  al.  (2008) 
that  govern  transmitted  fraction  in  a  canopy: 

Hc 

H*  —  T7~  (35) 

Mdc 


and 


T*  _  ^HldHcanl  _  H/ 

1  -  /k(Hc)~  /Hr  WW 

In  Eqs.  35  and  36,  Hc  is  the  height  of  the  canopy  and  Ibk  is  the  initial  height  of  the  dust  cloud  at 
the  roadside,  assumed  to  be  2  m  in  this  work.  A  height  of  2  m  was  selected  as  the  typical  height 
of  an  initial  vehicle  plume  based  on  the  previous  work  of  Veranth  et  al.  (2003).  The  sensitivity  of 
the  deposition  results  on  the  choice  of  H  is  shown  in  Fig.  7  in  Pardyjak  et  al.  (2008).  For  H* 
values  less  than  about  4,  the  results  are  quite  sensitive  to  this  choice,  however  the  simple  model 
described  here  allows  one  to  determine  a  range  of  potential  deposition  values  (based  on  a  range 
of  H*  values)  quickly  and  easily.  K(HC)  is  the  vertical  turbulent  diffusivity  at  the  top  of  the 
canopy  and  is  a  function  of  atmospheric  stability;  a  is  the  vegetative  surface  area  per  unit  volume, 
zt  and  Td  are  the  time  scales  of  vertical  turbulent  diffusion  and  of  particle  deposition  respectively. 
H*  detennines  IF  and  T  determines  the  canopy  removal  efficiency.  Here,  we  propose  to 
parameterize  IF  as  a  function  of  H 

IF  —  1  —  exp  (-bH*),  (37) 

where  b  is  a  positive  constant.  It  is  assumed  that  the  vegetative  canopy  is  located  at  the  roadside 
close  to  the  initial  plume  because  the  greatest  particle  filtering  should  occur  with  this 
configuration.  As  the  distance  is  increased  between  the  road  and  vegetation,  the  plume  can  mix 
vertically  and  the  effective  height  of  the  plume  compared  to  the  height  of  the  vegetation  will 
increase.  A  discussion  of  the  importance  of  the  distance  from  the  roadside  to  the  vegetation  can 
be  found  in  Pardyjak  et  al.  (2008)  and  Mao  et  al.  (2013).  The  form  of  Eq.  37  was  chosen 
considering  the  simulation  results  of  Pardyjak  et  al.  (2008).  Pardyjak  et  al.  (2008)  also  indicates 
that  the  canopy  removal  efficiency  approaches  1  as  7  becomes  large.  For  these  conditions  Eq. 
37  reduces  to  TF  =(1-IF).  Substituting  TF  =  c\p(-bll  )  into  the  equation  and  solving  for  IF 
yields  Eq.  37. 

The  canopy  removal  efficiency  depends  upon  T*.  To  calculate  T  parameterizations  for  a, 
Kz:(I I can),  and  Vd  are  needed.  y=LAI IHcan  where  LAI  is  the  leaf  area  index  and  is  physically 
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interpreted  as  the  vegetative  surface  area  per  ground  area.  For  many  locations,  LAI  can  be 
obtained  from  the  MODIS  satellite  database: 

(http://daac.oml.gov/cgibin/MODIS/GLBVIZ  1  Glb/modis  subset  order  global  col5.pl  last 

accessed  March  2013)  or  it  can  be  calculated  from  the  aerodynamic  roughness  according  to  Eq. 
4.6  (Choudhury  and  Monteith  1988) 

LAI  =  (zo/(0.28Hc))2/0.2.  (38) 

In  Eq.  38,  the  roughness  due  to  larger  scale  surface  irregularities,  such  as  sand  dunes,  is  assumed 
to  have  negligible  effect  upon  zQ  compared  to  the  vegetation.  To  calculate  Kzz(Hcan)  surface  layer 
similarity  is  assumed,  namely, 


Kzz(tfc)  =  KU,(ffc  -  d)/0«7/c  -  d)/L), 


(39) 


where  k  is  the  von  Karman  constant  =0.4,  L  is  the  Monin-Obukhov  length  scale  defined  as 
L  =  u» /\rcQa(g  /  T0^j\.  To  and  Q0  are  the  absolute  surface  temperature  and  kinematic  heat  fluxes 

respectively,  u*  (the  friction  velocity)  is  the  turbulent  velocity  scale,  and  (j)(H-d)/L  is  the 
universal  stability  function  are  given  by  Arya  (2001)  as: 
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where  d  is  the  displacement  height  of  the  canopy.  It  can  be  detennined  from  the  vertical  profile 
of  u  (the  temporally  averaged  velocity)  or  modeled  (Pardyjak  et  al.  2008).  The  deposition 
velocity,  Vd,  is  a  parameter  to  quantify  removal  of  particles  by  impaction  (Seinfeld  and  Pandis 
1998).  It  has  been  measured  experimentally  by  Zhu  et  al.  (2011)  and  is,  in  general,  a  function  of 
the  modes  of  particle  removal:  sedimentation  (due  to  gravitational  settling),  interception  (a  result 
of  when  flow  streamlines  transport  a  particle  into  contact  with  a  vegetative  surface),  Brownian 
motion,  and  impaction  (a  result  of  when  the  time  scale  of  the  flow  disturbances  resulting  from  a 
vegetative  element  is  much  smaller  than  the  time  scale  of  the  particle’s  response  to  the  changing 
flow)  (Seinfield  and  Pandis  1998).  It  is  difficult  to  parameterize  Vd  ;  many  impaction  models 
require  estimates  for  the  size  and  orientations  of  vegetative  elements  (e.g.  Aylor  and  Flesh  2001; 
Slinn  1982),  which  are  unknown  for  the  three  field  studies  and  are  difficult  to  parameterize. 
Because  of  the  difficulty  of  estimating  Vd  ,  an  approximation  for  Vd  is  utilized  in  this  work  to 
parameterize  T  .  It  is  assumed  that  Vd  is  proportional  to  the  turbulent  velocity  scale  u*.  The  basis 
for  this  assumption  is  found  by  examining  the  Stokes  number,  Stk.  The  Stokes  number  describes 
the  ability  of  a  particle  to  faithfully  follow  a  streamline  in  a  flow,  and  may  be  thought  of  as  a 
dimensionless  time  scale  that  describes  the  likelihood  of  a  particle  to  deposit  onto  an  impaction 
surface  (Hinds  1982).  A  common  parameterization  used  in  field  studies  is  (Seinfeld  and  Pandis 
1998): 
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Stk  =  (u*vs)/(gR), 


(42) 


where  vs  is  the  gravitational  settling,  g  is  the  acceleration  due  to  gravity,  and  R  is  the  collector 
radius.  Assuming  that  Vd  is  proportional  to  Stk  and  that  vs  and  R  are  constant  for  the  three  field 
studies,  the  conclusion  that  Vd  is  directly  proportional  to  u*  results.  This  is  expressed  as  the  new 
parameter  T ^ 


(43) 


Because  field  measurements  of  Vd  indicate  that  it  is  smaller  than  u*  (Zhu  et  al.  2011)  it  is 
assumed  that  7^  =£  T*.  However,  for  canopies  of  identical  mean  vegetative  element  size  and 

•  ••••a  .  %  f 

PMio  plumes  of  identical  size  distributions  Tm  and  T  are  proportional  to  each  other. 

Tm  can  be  calculated  based  on  field  or  model  data.  For  situations  in  which  neither  is  available, 
the  following  methodology  can  be  used: 

1.  Use  Table  10.5  in  Arya  (2001)  to  determine  z0  for  the  land  type  (grass,  forest). 

2.  Use  the  equation  y  =  (zo/(0.28/f.1,5)2)/0.2  to  calculate  y  (the  vegetative  surface  area 
per  unit  volume  of  vegetation. 

3.  Use  a  meteorological  network  such  as  mesowest  (url  http://mesowest.utah.edu)  to  obtain 
the  mean  wind  speed,  u  at  a  reference  height  of  10  m. 

4.  Use  Seinfield  and  Pandis  (2001)  to  obtain  the  Monion-Obukhov  length  scale,  L  for 
cloudy  conditions  L. 

5.  Use  Eq.  39  to  calculate  u*//f(//c). 

This  methodology  can  be  applied  by  practitioners  to  detennine  characteristics  of  the  vegetation 
needed  to  achieve  a  desired  PMio  reduction  in  the  nearest  100  m  downwind  of  an  unpaved  road. 
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4.  Results  and  Discussion 


4.1.  Laboratory  Wind  Tunnel  Measurements 

4.1.1.  Deposition  Experiments 

The  Stokes  number  (Stk)  is  a  dimensionless  parameter  that  describes  the  likelihood  of  a  particle 
to  deposit  onto  an  impaction  surface  (Hinds  1982).  This  can  also  be  understood  as  the  deposition 
efficiency.  Deposition  efficiency  is  the  ratio  of  deposited  particles  onto  an  impaction  surface 
over  the  total  volume  of  air  occupied  by  the  frontal  impaction  surface  area  per  unit  time.  The 
Stokes  number  is  shown  below  in  Eq.  44,  where  pp  is  particle  density,  dp  is  particle  diameter,  Cc 
is  the  Cunningham  correction  factor  (~1 ),  pair  is  the  dynamic  viscosity  of  air,  and  Ls  is  the 
effective  length  scale  of  the  impaction  surface. 


Stk  = 


p  cfCU 

*  p  p  c  o 

18//  L 

t  air  s 


(44) 


In  order  to  cover  a  larger  Stk  range  and  acquire  data  from  the  limit  of  detection  (0%)  to  the  point 
of  saturation  (100%),  three  substrate  sizes  and  three  velocities  were  used  (i.e.,  Ls  =1.4,  1.0,  and 
0.5  cm  and  Ua  =  2,  5,  and  8  m/s).  Regardless  of  the  variation  in  these  parameters,  our  data  for  DF 
versus  Stk  appear  to  be  almost  identical  to  the  shape  of  the  curve  developed  from  the  experiments 
of  Hinds  (1982).  The  difference  is  that  our  data  have  much  smaller  Stk  (due  to  relatively  small 
particles)  and  our  flow  was  turbulent.  However,  our  present  data  indicate  significant  deposition, 
ranging  from  0%  to  100%,  but  at  much  smaller  Stk  than  expected  based  on  Hinds  (1982).  This 
seems  to  suggest  our  hypothesis  that  turbulent  mixing  enhances  deposition  onto  surfaces  is  valid. 
Furthermore,  R/1  can  be  used  as  the  appropriate  parameter  to  scale  the  Stk  and  fit  to  the  curve 
developed  from  Hinds. 

Figure  25  shows  (i)  a  plot  of  DF  versus  Stk  for  the  classical  laminar  flow  experiments  from 
Hinds  (1982)  for  particle  impaction  onto  a  flat  plate,  (ii)  laminar  flow  data  from  the  experiments 
of  May  and  Clifford  (1967)  using  a  vertical  wind  tunnel  (neglecting  gravitational  settling),  and 
(iii)  grid-generated  turbulence  data  from  our  present  experiments.  Also  shown  in  this  plot  is  our 
turbulence  data  that  were  scaled  by  R^0  4.  The  reasoning  for  scaling  by  the  Taylor-scale 
Reynolds  number  has  been  described  in  the  previous  section. 
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Figure  25:  Plot  of  DF  (deposition  fraction)  as  a  function  of  Stokes  numbers  where  the  solid  black  line  is 
the  Hinds  data  (1982),  the  green  dashed  line  is  from  the  May  and  Clifford  experiments  (1967),  (red 
circles)  are  from  the  present  experiments,  and  (blue  diamonds)  are  our  present  experiments  scaled  by 

R,0'34. 

Figure  26  illustrates  the  decay  of  turbulence  measured  in  our  wind  tunnel  compared  to  a  similar 
grid-turbulence  experiment  by  Comte-Bellot  and  Corrsin  (1965).  The  turbulence  intensity  (ut)  at 
each  point  is  measured  from  the  time  series  hot-wire  data  using  Eq.  10: 


u . 

In  Eq.  45,  u  is  the  mean  velocity  component  in  the  streamwise  direction  and  u  ’  is  the  r.m.s.  of  the 
mean  velocity.  In  theory,  our  turbulence  intensities  for  every  velocity  measurement  should  scale 
exactly  with  each  other  at  every  point  measured  due  to  the  behavior  of  isotropic  turbulence. 
However,  there  appears  to  be  variation  in  the  calibration  curves  from  before  and  after  the  data 
collection.  Review  of  the  two  calibration  curves  indicated  that  there  was  a  slight  variation 
between  velocity  measurements  before  and  after  data  collection.  The  calibration  curves  are  third 
order  polynomial  fits  of  the  voltage  from  the  hot-wire  to  the  velocity  derived  from  the  pressure 
measurements  under  the  idealized  conditions  in  calibration  unit.  In  order  to  account  for  the 
variation  between  the  curves,  the  two  polynomial  fits  were  averaged  and  a  single  equation  was 
used  for  data  processing  (R2  ~  0.998).  The  use  of  this  average  fitted  equation  would  suggest 
errors  might  have  propagated  in  the  calculations  of  the  instantaneous  velocity  and  further  into  the 


44 


turbulence  intensity  calculations.  Figure  26  demonstrates  the  decay  of  turbulence  intensities  for 
our  grid  built  in  the  University  of  Utah  Environmental  Fluids  Lab.  The  different  experiments 
from  the  Comte-Bellot  and  Corrsin  (1965)  and  our  present  data  follow  slightly  different  slopes, 
owing  to  the  different  grid  solidities  used  (0.44  for  Comte-Bellot  and  Corrsin  and  0.50  for  our 
present  data). 

Note  in  Figure  26,  that  the  distance  from  the  grid  (x),  is  nonnalized  by  the  mesh  spacing  of  the 
turbulence  grid  ( M ). 


Figure  26:  Plot  showing  the  decay  of  normalized  grid-generated  turbulence  in  the  University  of  Utah 
wind  tunnel.  Present  data  follows  a  line  proportional  to  (x/M)'0'9  and  the  Comte-Bellot  and  Corrsin  data 
follows  a  line  proportional  to  (x/M)'1'3. 


Considering  DF  results,  the  same  DF  curve  at  a  lower  Stokes  number  indicates  that  turbulent 
motion  appears  to  increase  the  DF  onto  the  impaction  surface.  However,  this  does  not  indicate  if 
the  turbulent  motions  are  enhancing  deposition  on  the  other  axes  as  well.  Figure  27  illustrates 
the  effect  of  isotropic  turbulence  on  deposition  for  all  six-deposition  surfaces.  Clearly,  grid 
turbulence  results  in  a  significant  deposition  increase  on  all  axes.  It  should  be  noted  that  the  ratio 
of  DF  for  turbulent  to  laminar  flow  (i.e.  turbulence  grid  to  no  turbulence  grid)  is  considerably 
more  on  the  non-impaction  surfaces.  This  appears  to  indicate  that  turbulent  mixing  has  a  greater 
influence  on  the  non-impaction  surfaces  (Moran  et  al.  2011).  Further,  deposition  onto  the  surface 
that  is  aligned  with  gravitational  vector  (+z,  where  gravitational  settling  would  be  important)  is 
not  significantly  different  from  the  other  surfaces.  This  counter  to  what  would  be  expected  for 
larger  particles  (see  e.g.,  Seinfeld  and  Pandis,  1998),  where  deposition  onto  the  +z  surface  would 
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be  very  important.  Comparing  the  overall  deposition  between  grid  and  no-grid  experiments,  there 
was  a  300%  increase  in  deposition  when  the  grid  was  used. 


Laminar  Flow  without  Grid 
■  Grid  Generated  Turbulence 


1 1  i  1 1 

-x  +x  -y  +y  -z  +z 

Axis  of  Deposition 

Figure  27:  Aerosol  DF  for  each  of  the  various  substrate  axes  for  laminar  and  turbulent  flow  (U„  =  5  m/s 
and  Ls  —  1  cm). 

The  previous  results  have  been  utilized  to  parameterize  deposition  fraction,  DF,  for  turbulent 
conditions.  While  doing  this  a  new  parameter  is  introduced,  a  modified  Stokes  number,  Stk*.  It  is 
the  product  of  the  traditional  stokes  number  and  the  Taylor-scale  Reynolds  number,  R,,  This  is 
expressed  as: 

Stk*  =  StkR%.  (46) 

Stk  is  utilized  because  of  the  results  of  the  wind  tunnel  measurements  suggest  that  DF  is 
negligible  for  Stk*  <  1  and  becomes  significant  as  Stk  increases  from  1 .  The  traditional  Stk  does 
not  yield  a  simple  conclusion.  DF  is  related  to  Stk  by  the  following  expression: 

DF=(  1 00- 1 00/(440. 5  *Stk*)3M)  (47) 

Eq.  47  has  important  impacts  when  using  deposition  instead  of  the  traditional 
D  F=(Stk/( S tk+0 . 8 ) )2 .  Both  expressions  require  an  estimate  of  mean  vegetation  element  size  to 
estimate  DF.  The  traditional  expression  requires  the  mean  vegetative  element  size,  de,  to  be 
extremely  small,  on  the  order  of  50  pm  to  give  good  agreement  with  PMio  field  data.  In  contrast 
Eq.  47  requires  a  more  reasonable  5-15  mm.  Eq.  47  is  also  less  sensitive  to  small  changes  to  de. 
The  collapsed  results  with  the  new  parameterization  are  show  in  Fig.  28  below. 
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Amatul  (2006)  implemented  the  traditional  deposition  model  into  QUIC  for  canopies.  One 
important  concern  with  her  model  was  the  need  for  very  small  vegetative  length  scales  to 
produce  experimentally  observed  deposition  levels.  This  current  work  has  replaced  the  Amatul 
(2006)  model,  which  required  unrealistically  small  values  of  de  to  obtain  plausible  results  with 
the  new  model,  which  has  realistic  values  of  de  that  correspond  to  field  study  results  (See  Section 
3.2). 


Figure  28:  Deposition  fraction,  DF,  as  a  function  of  Stk*  showing  the  collapse  of  the  artificial  substrate 
data. 

4.1.2.  Vegetation  Deposition  Experiment 

In  addition  to  the  above  artificial  substrate  deposition  experiments,  artificial  vegetation,  twig 
elements,  and  rods  were  tested  to  verify  that  the  Stk  scaling  was  valid  for  other  geometry.  A 
total  of  60  wind  tunnel  tests  were  performed  (after  quality  controls)  on  7  size  categories  and 
wind  speeds  of  4.9,  8.3,  and  9.5  m/s.  Generally  the  results  agree  with  the  Stk  parameterization, 
however,  it  is  clear  that  specific  geometry  and  the  methods  to  calculate  deposition  fraction  affect 
the  parameterization. 

The  vegetation  deposition  fraction  results  clearly  indicate  that  the  classical  laminar  theory  Stokes 
number  underestimates  the  dry  deposition  for  all  vegetation  and  rod  elements.  This  result  is  the 
same  as  the  artificial  substrates  also  shown  in  Figure  29. 
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Figure  29:  Deposition  fraction,  DF,  as  a  function  of  the  classical  laminar  Stokes  number  for  all 
experiments  including  vegetation.  Deposition  values  are  shown  for  broad  leaf  (flat  into  the  wind),  edge 
leaf  (broad  leaf  edge  onto  the  wind),  artificial  onion  grass,  twigs,  and  aluminum  and  brass  welding  rod  of 
several  diameters.  Error  bars  indicate  the  standard  deviation  for  identical  experiments.  The  artificial 
substrate  results  and  parameterization  are  shown  in  red. 

Applying  the  modified  Stokes  number  described  by  Eq.  46,  it  appears  that  the  parameterization 
fit  is  good  for  low  Stk  values  and  loosens  for  higher  Slk  as  shown  in  Figure  30.  This  may  be 
related  to  the  reason  the  deposition  fraction  appears  greater  than  100%.  Eq.  1  is  used  to  calculate 
the  deposition  fraction  based  on  the  air  velocity,  the  aerosol  concentration,  the  length  of  time, 
and  the  frontal  area  of  the  surface.  Basically,  the  deposition  fraction  is  the  mass  deposited  on  the 
surface  over  the  mass  of  particles  in  the  volume  of  air  that  the  surface  has  seen  projected  through 
time.  In  the  case  of  the  vegetation  tests,  DF  values  in  excess  of  100%  are  possible  if  the  Stokes 
number  is  high  enough  to  cause  very  high  impaction  on  the  windward  surface.  Any  additional 
turbulence  enhanced  deposition  on  the  downwind  surface  will  be  added  to  the  deposited  mass, 
thus  possibly  producing  DF  >  100%.  The  results  showing  DF  greater  than  100%  serve  to 
illustrate  the  effects  of  downwind  deposition  and  to  illustrate  the  deficiency  of  using  the 
projected  area  when  calculating  deposition  fraction. 
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Figure  30:  Deposition  fraction,  DF,  as  a  function  of  the  modified  Stokes  number  (Stk*).  Error  bars 
indicate  the  standard  deviation  for  identical  experiments.  The  artificial  substrate  results  and 
parameterization  are  shown  in  red. 

The  results  show  that  the  broad  leaf  elements,  both  faced  into  the  wind  and  oriented  edge  into  the 
wind,  achieve  low  deposition  rates  as  the  Stokes  number  is  very  low.  The  other  elements 
including  grass  and  rods  show  much  higher  deposition  as  expected  by  the  higher  Stokes  number 
(and  Stk  ). 

Several  alternative  scaling  parameters  were  investigated  to  produce  a  better  fit,  but  none  were 

found  that  produced  satisfactory  results.  It  is  believed  that  there  are  a  combination  of  factors  that 

may  be  significant  including  the  ratio  of  the  Taylor  Microscale  to  the  characteristic  length  (A/D) 

•  •  •  •  ■  ^ 
and  the  geometry  specific  width  distributions  where  one  element  has  many  Stk  values.  In  order 

to  keep  the  parameterization  as  general  as  possible  (not  species  specific),  the  authors  have 

chosen  to  incorporate  the  parameterization  from  the  artificial  substrates  (Eq.  47)  into  the  QUIC 

model  (Section  3.2). 
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4.2.  Comparison  of  the  Improved  QUIC  model  to  the  Full-scale  Deposition  Field 
Experiment  Results 


4.2.1.  Hanford,  Washington  Semi-arid  Vegetation  Experiment 

In  this  section  results  from  the  Hanford  site  field  experiments  are  compared  with  the  improved 
QUIC  model.  For  these  simulations,  the  coordinate  system  was  defined  as  shown  in  Figure  1. 
That  is,  Y  was  parallel  to  road,  X  was  perpendicular  to  the  road,  and  Z  was  vertical.  The 
computational  domain  was  1000  m  x  180  m  x  50  m  in  the  Y,  X  and  Z  directions  respectively.  A 
fetch  upwind  of  the  road  was  specified  to  be  30  m.  The  computational  grid  was  uniform  in  each 
direction  and  given  by:  dy  =  200  m,  dx  =  1  m,  and  dz  =  0.25  m.  The  vegetation  height  in  QUIC 
is  uniform  and  specified  to  be  1 .4  m,  which  was  the  average  vegetation  height  measured  in  the 
field.  For  all  simulations  1,000,000  particles  of  5  pm  diameter  were  released.  QUIC  also  requires 
an  attenuation  coefficient  (e.g.  Cionco  1965)  to  model  the  velocity  profile  within  the  vegetative 
canopy  as  described  by  Pardyjak  et  al.  (2008).  The  attenuation  coefficient  was  selected  to  best 
match  the  observed  velocity  profile  (see  Fig.  32).  Within  the  vegetation,  an  attenuation 
coefficient  of  1 .4  was  specified.  The  input  meteorological  conditions  for  each  QUIC  simulation 
were  specified  from  the  experimental  data  presented  in  Table  2. 

Objective  1 

The  atmospheric  conditions  for  the  June  7,  2011  Objective  1  IOP  at  Handford,  WA  were 
unstable  as  indicated  by  a  Monin-Obukhov  length  scale  L  =  -47.7  m.  The  average  turbulent 
kinetic  energy  was  0.4  m2/s2  and  the  turbulent  velocity  scale  a*  was  0.61  m/s.  All  of  these 
quantities  were  evaluated  from  measurements  made  at  a  height  of  3.9  m  (which  should  be  well 
within  the  constant  flux  layer  of  the  atmosphere).  The  aerodynamic  roughness  height  zQ  =  0.07  m 
and  the  average  wind  speed  at  3.9  m  was  5.4  m/s  during  the  experiment. 

The  experimentally  measured  horizontal  flux,  along  with  its  location,  and  the  calculated  Fmin 
value  utilizing  Eq.  23  are  shown  in  Table  7  and  Figure  31  where  the  horizontal  flux  F  is 
nonnalized  by  the  roadside  flux,  Froadside. 


Tower 

Location 
Downwind  (m) 

Measured 

Flux 

F/Froad  sj(/e 

Measured 

F mit/F roadside 

QUIC 

Estimated 

Fmii/Froad  side 

Tower  1 
(Roadside) 

0 

1.00 

0.72 

0.71 

Tower  2 
(Middle) 

6.4 

0.86±0.14 

0.72±0.1 

0.71 

Tower  3 
(Downwind) 

56.4 

0.74±0.13 

0.71±0.1 1 

0.71 

Table  7:  Results  from  the  Objective  1  IOP  on  June  7,  2011  (unstable  atmospheric  conditions).  QUIC 
simulations  predicted  an  Fmir/Froadside  value  of  0.71  given  the  site  roughness  and  meteorological  conditions. 
The  measured  values  in  the  right  most  column,  closely  agree  with  the  QUIC  value. 
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Figure  31:  The  measured  decrease  in  normalized  horizontal  PMi0  flux,  F,  compared  to  the  QUIC 
simulation  of  the  same  of  the  Handford,  WA  June  7,  201 1  IOP.  Bars  are  95  %  confidence  intervals. 

Again,  the  first  objective  of  the  flux  study  was  to  validate  the  QUIC  deposition  model  for 
vegetative  canopies.  The  quantity  Fmin/Froadside  was  utilized  as  a  metric  to  estimate  the  accuracy 
of  QUIC  (See  Fig.  3 1).  The  QUIC  computations  produced  a  value  of  0.7 1  for  F min/F roads\de-  This 
may  be  compared  to  the  measured  results  in  the  right  most  column  of  Table  7.  The  measured 
results  were  calculated  using  Eq.  22  with  the  middle  and  far  downwind  towers.  The  result  was 
then  normalizing  the  by  the  measured  Froadside-  The  relative  difference  between  the  middle  towers 
measured  0.72  and  the  QUIC  estimated  0.71  is  only  1.4%  suggesting  that  for  the  June  7  IOP 
QUIC  is  modeling  the  processes  of  deposition  and  dispersion  well. 

Objective  2 

The  atmospheric  conditions  for  the  June  13,  2011  Objective  2  IOP  at  Handford,  WA  were 
weakly  unstable  as  indicated  by  a  Monin-Obukhov  length  scale  L  =  -1 17  m.  In  the  constant  flux 
layer  above  the  canopy,  the  average  turbulent  kinetic  energy  was  0.3  m  /s  and  the  turbulent 
velocity  scale  w*  was  0.71  m/s.  All  of  these  quantities  were  evaluated  from  measurements  made 
at  a  height  of  3.9  m  (which  should  be  well  within  the  constant  flux  layer  of  the  atmosphere).  The 
aerodynamic  roughness  height  z0  =  0.07  m  and  the  average  wind  speed  at  3.9  m  was  6.8  m/s 
during  the  experiment.  The  QUIC  models  for  mean  advection  and  turbulent  kinematic  shear 

stress  it*  =  V it'w'  within  canopies  are  shown  in  Figures  32  and  33.  The  model  for  the  mean 
wind  shows  excellent  agreement.  The  model  for  turbulent  shear  stress  is  satisfactory;  however 
the  QUIC  model  displays  a  non-physical  corner  at  the  top  of  the  canopy.  This  is  most  likely  a 
consequence  of  the  simple  canopy  model  that  assumes  a  constant  vegetative  density,  y,  with 
height  within  the  canopy. 
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Figure  32:  Comparison  of  the  QUIC  vegetative  canopy  mean  wind  model  to  3D  sonic  anemometer  data 
taken  37  nr  downwind  of  the  road,  during  the  June  13,  201 1  10P  at  Hanford,  WA. 


Figure  33:  Comparison  of  the  QUIC  canopy  turbulence  model  to  3D  sonic  anemometer  data  37  m 
downwind  of  the  road,  taken  during  the  June  13,  201 1  IOP  at  Handford,  WA.  Here  u*  is  a  local  turbulent 
kinematic  shear  stress. 


52 


- QUIC 

■  Handford  Data 


Figure  34:  The  variation  of  concentration  with  height  within  the  canopy  at  Hanford,  WA.  The  field  data 
indicates  that  concentration,  C  increases  with  height  in  an  irregular  manner  that  is  generally  followed  by 
the  QUIC  simulation  except  near  the  ground.  These  data  were  taken  32  m  downwind  of  the  road.  Error 
bars  indicate  95%  confidence  intervals. 

Figure  34  shows  the  high-resolution  concentration  profile  taken  32  m  downwind  of  the  road. 
The  simulation  and  field  data  show  good  agreement  except  near  the  ground.  This  is  most  likely  a 
consequence  of  deposition  to  the  ground  being  underestimated.  The  QUIC  simulation  shows  a 
decrease  in  concentration  with  increasing  height  in  the  bottom  half  of  the  canopy.  In  the  top  half 
of  the  canopy  QUIC  predicts  an  increase  in  concentration  with  increasing  height.  The  data  show 
the  same  basic  trends  but  discrepancies  between  simulation  and  data  occur  at  multiple  points. 
This  is  most  likely  due  to  canopy  heterogeneities.  To  better  understand  the  significance  in  the 
variability  shown  in  Fig.  34,  we  note  that  the  colocation  tests  had  a  standard  deviation  of 
approximately  0.06  mg/m3  (See  Appendix  A).  In  Fig.  34,  the  normalizing  concentration  was 
C|z=i.3Hc  =  41  mg/m  .  Hence,  the  differences  in  Fig.  34  which  are  greater  than  -0.06/41  =  0.0015 
are  larger  than  the  inter-instrument  variability.  This  indicates  the  curvature  in  the  measurements 
could  indeed  be  real. 
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To  better  quantify  the  predicted  error  from  the  QUIC  simulations,  we  introduce  the  following 
metrics  that  are  often  used  (Warner  et  ah,  2004)  to  evaluate  dispersion  models.  In  the  following 
equations  Cp  indicate  predicted  values  from  running  QUIC,  while  C0  are  observed  measurements 
from  the  field  experiments.  The  following  metrics  have  been  computed: 

FB  -  Fractional  Bias 


0.5(Cp  +  C0) 

RE  -  Relative  Error 

RE  =  Cp_C°  *  100% 
C0 


NMSE  -  Normalized  Mean  Square  Error 

NMSE  = 

Cp  C0 

Here,  the  overbar  indicates  an  average  over  all  of  the  samples.  The  first  two  metrics  are 
indicative  of  the  model’s  ability  to  over  predict  or  under  predict  a  value,  while  the  NMSE  is 
measure  of  the  spread.  Figure  35  Below  illustrates  these  paired  comparison  parameters,  while 
Table  8  shows  the  numeric  values  of  the  parameters.  Here,  Mean  Wind  is  a  comparison  of  mean 
wind  speed,  Turbulence  is  a  comparison  of  the  local  u*  and  Concentration  is  local  concentration. 


Quantity 

FB 

NMSE 

RE  (%) 

Mean  Wind 

-0.00063 

0.00621 

0.24999 

Turbulence 

0.13993 

0.02955 

15.639 

Concentration 

0.10187 

0.02687 

11.682 

Table  8:  Error  metrics  for  the  Handford,  WA  experiment  for  the  mean  wind,  turbulence  and 
concentration  measurements. 
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Figure  35:  Normalize  scatter  plots  comparing  the  QU1C  predicted  mean  wind  speed,  turbulence  and 
concentrations  with  the  field  observations  taken  at  Handford,  WA. 
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The  heterogeneity  of  the  canopy  is  shown  in  the  results  from  the  LAI  measurements  in  Figure  36. 
The  measurements  indicate  that  the  majority  of  the  vegetative  surface  area  is  located  near  the 
ground.  For  example  the  top  half  of  the  canopy  contains  only  about  10%  of  the  vegetative 
surface  area.  This  is  in  contrast  to  the  simple  canopy  model  that  has  been  implemented  in  QUIC, 
which  assumes  that  vegetative  surface  area  is  distributed  evenly  vertically  and  horizontally.  This 
plot  indicates  there  is  tremendous  variation  of  vegetative  density  in  both  the  vertical  and 
horizontal.  Despite  heterogeneities  present  in  the  canopy,  the  simple  homogenous  canopy  model 
offers  a  simple  modeling  tool  for  estimating  the  transport  and  removal  of  PMio  in  conditions 
typical  of  military  ranges.  The  error  metrics  presented  above  provide  an  initial  indication  of  the 
uncertainty  in  the  model.  Further  experimental  data  are  necessary  to  understand  the  models 
ability  to  be  widely  applicable  over  different  types  of  vegetation. 
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LAI/LAIground 

Figure  36:  The  decrease  of  LAI  with  height  within  the  canopy  at  Hanford.  The  figure  indicates  the 
majority  of  the  vegetative  surface  area  is  located  near  the  ground.  For  example  the  plot  indicates  100%  of 
the  vegetative  surface  area  is  located  above  z/Hc  =  0,  slightly  less  than  40%  is  located  above  z/Hc  =  0.25, 
and  so  forth.  Also  note  the  tremendous  variation  in  LAI  measured  as  noted  by  the  standard  deviation  error 
bars.  These  data  were  taken  from  48  LAI  profiles  taken  on  a  uniform  grid  with  the  3D  sonic  tower  being 
the  origin. 
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Impact  of  the  New  Deposition  Model 

Figure  37  shows  the  dependence  of  the  simulated  QUIC  transmitted  fractions  ( TF  =  Fmin/Froad 
side)  on  the  vegetative  element  size  that  must  be  specified  in  the  simulations,  de,  for  the  Handford, 
WA  site  on  June  7,  201 1.  Figure  37a  shows  a  traditional  removal  model  that  neglects  turbulence 
enhancement  of  PMio  deposition  and,  while  Fig.  37b  shows  results  from  the  new  model  (Eq.  47) 
that  incorporates  turbulence  enhancement.  Note,  that  while  TFs  with  reasonable  values  can  be 
achieved  with  standard  models  in  which  turbulence  enhancement  is  neglected,  unrealistically 
small  vegetation  elements  (~50  pm)  must  be  used  to  explain  the  measured  removal  rates 
(TF~70%)  measured.  If  turbulence  effects  are,  included  realistic  values  for  (dt,  ~  10  mm)  explain 
the  removal  measured.  In  summary  vegetation  in  canopies  (and  windbreaks)  have  higher  PMio 
removal  rates  than  suggested  by  traditional  removal  theory.  The  Flandford  June  7,  2011  data 
indicate  that  the  new  particle  deposition  model  introduced  into  QUIC  as  part  of  this  project  has 
greater  physical  fidelity  than  traditional  models  such  as  those  used  by  Amatul  (2006)  that  utilize 
Eqs.  16  and  17.  Note  also  that  the  new  model  also  shows  smaller  sensitivity  to  de  as  well.  This 
has  large  implications  given  the  relatively  high  uncertainty  which  exists  in  estimating  vegetative 
element  sizes. 


a)  b) 


Figure  37:  Dependence  of  QUIC  measured  transmitted  fraction  upon  the  vegitative  element  size,  de, 
assumed  in  (a)  the  orginal  QUIC  model  and  (b)  the  modified  QUIC  model  shown  for  the  Handford,  WA 
June  7,  2011  field  study.  Panel  a)  displays  results  from  a  traditional  PM10  removal  model,  that  neglects 
turbulence  enhancement  to  deposition  (Amatul  2006).  In  panel  b)  the  new  deposition  model  developed  in 
this  work  includes  turbulent  enhancement  of  deposition.  Accounting  for  the  effects  of  turbulence  on 
deposition  results  in  realistic  estimates  for  de ;  neglecting  it  requires  unrealistically  small  values  for  de  to 
account  for  measured  removal.  Transmittance  Fraction  is  defined  Section  2.5. 


57 


4.2.2.  Raft  River  Windbreak  Experiment  -  Malta,  Idaho 

Because  of  two  equipment  failures,  the  wind  data  and  PMio  concentration  data  were  taken  during 
two  separate  IOPs.  The  wind  data  were  taken  during  the  July  12th  IOP  and  the  concentration  data 
were  taken  on  the  July  27th  IOP.  Even  though  a  failure  occurred  that  deleted  the  majority  of  the 
3D  sonic  wind  data  on  the  July  27th  IOP,  the  upstream  3D  sonic  data  survived.  This  is  important 
because  these  data  are  necessary  to  generate  an  accurate  QUIC  simulation  of  the  PMio  dispersion. 

The  two  thick  windbreaks  each  had  an  optical  porosity  of  approximately  0.01  and  fetches 
(thicknesses)  of  7  m.  The  optical  porosity  was  determined  by  first  computing  an  aerodynamic 
porosity  (a)  from  the  3D  sonic  anemometer  data  subsequently  using  the  a=P°  4  model  to  solve  for 
P,  the  optical  porosity.  The  height  of  each  windbreak  was  approximately  4  m.  The  upstream 
roughness  was  0.02  m. 

In  this  section  results  from  the  Raft  River  Windbreak  field  experiments  are  compared  with  the 
improved  QUIC  model.  For  these  simulations,  the  coordinate  system  was  defined  as  shown  in 
Figure  1.  That  is,  Y  was  parallel  to  road,  X  was  perpendicular  to  the  road,  and  Z  was  vertical. 
The  computational  domain  was  1000  m  x  180  m  x  50  m  in  the  Y,  X  and  Z  directions  respectively. 
A  fetch  upwind  of  the  road  was  specified  to  be  30  m.  The  computational  grid  was  uniform  in 
each  direction  and  given  by:  dy  =  200  m,  dx  =  1  m,  and  dz  =  0.25  m.  For  all  simulations 
1,000,000  particles  of  5  pm  diameter  were  released  were  released.  The  input  meteorological 
conditions  for  each  QUIC  simulation  were  specified  from  the  experimental  data  presented  in 
Table  2. 

July  12th  IOP 

This  nighttime  IOP  measured  the  mean  and  turbulent  winds  created  by  the  presence  of  the 
windbreak.  The  atmospheric  conditions  were  slightly  stable  with  a  Monin-Obukhov  length  scale, 
F  =  95  m  that  was  measured  at  a  height  of  3.25  m  upstream  of  the  windbreaks.  The  turbulent 
velocity  scale,  u*,  was  0.33  m  s'1  again  measured  at  3.25  m  upstream.  The  wind  speed  data  for 
winds  that  were  approximately  perpendicular  to  the  windbreaks  are  shown  for  Towers:  B  (x  = 
13.8  m),  C(x=24.8  m),  D(x=3.5  m),  and  E(x=55.8)  are  shown  in  Figures  38  and  39  along  with 
the  QUIC  simulation  results. 
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x  =  13.  8  m  (in  between  windbreaks) 


x  =  27.  8  m  (leeward  of  2nd  windbreak) 


x  =  42.1m  (3.5  HWb  Downwind  of  2nd  windbreak)  x=  55.4  m  (6.  75  HWb  Downwind  of  2nd  windbreak) 


Figure  38:  Comparison  of  QUIC  simulations  and  Raft  River  Windbreak  field  experiment  data  for  the 
mean  winds  from  the  July  12,  2012  10P.  Reference  values  measured  at  a  height  of  3.25  m  upstream  of 
wind  breaks. 
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x  =  13.  8  m  (in  between  windbreaks) 


x  =  27.8m  (leeward  of  2nd  windbreak) 


Figure  39:  Comparison  of  QUIC  simulations  of  windbreak  turbulence  (local  kinematic  turbulent  shear 
stress)  to  Malta  wind  10P  (July  12,  2012).  Reference  values  measured  at  a  height  of  3.25  m  and  35  m 
upstream  of  the  windbreaks. 

The  figures  indicate  that  the  QUIC  simulation  and  the  Raft  River  Windbreak  wind  data  is  able  to 
follow  the  trends  of  both  the  mean  wind  and  the  turbulence.  The  quiet  zone  is  evident  between 
the  windbreaks  (upper  left  panel)  as  well  as  the  high  shear  located  leeward  of  the  second 
windbreak  (upper  right  panel).  The  initial  recovery  of  the  flow  from  the  disruption  of  the 
windbreak  is  evident  in  the  lower  two  panels.  The  sensor  array,  which  was  limited  in  size  by  a 
highway,  was  not  sufficiently  large  to  capture  the  complete  recovery  of  the  flow  from  the 
disruption  of  the  windbreak. 

QUIC  suggests  some  trends  that  were  not  observed  by  our  sensor  array.  Namely,  concerning  the 
mean  wind,  the  lower  right  panel  of  Fig.  39  shows  a  slightly  elevated  velocity  “bump”  near  the 
ground  in  the  vertical  profile  of  u.  This  is  a  result  of  the  high  shear  zone,  which  starts  at  the  top 
of  the  windbreak,  reaching  the  ground  and  accelerating  the  local  flow.  Such  phenomena  are 
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difficult  to  observe  in  field  studies.  They  have,  however,  been  observed  in  wind  tunnel  studies 
such  as  Judd  et  al.  (1996).  Evidently,  the  new  QUIC  model  is  able  to  resolve  such  features. 

July  27th  IOP 

This  daytime  IOP  measured  dispersion  of  PMio  with  the  flow  the  presence  of  the  windbreak. 
The  atmospheric  conditions  were  very  unstable  with  a  Monin-Obukhov  length  scale  of  L  =  -2  m 
that  was  measured  at  a  height  of  3.25  m  upstream  of  the  windbreaks.  The  turbulent  velocity  scale, 
14,  was  0.36  m  s'1  again  measured  at  3.25  m  upstream.  The  measurements  at  Towers:  B  (x  =  13.8 
m),  C  (x=24.8  m)  and  D  (x  =  3.5  m),  are  shown  in  Figure  40. 


x  =  13.  8  m  (in  between  windbreaks)  x  =  27.  8  m  (leeward  of  2nd  windbreak) 


Figure  40:  Comparison  of  QUIC  simulation  results  to  the  Raft  River  Windbreak  (Malta,  ID) 
PMio  concentrations  for  the  July  27,  2012  IOP.  CS0urceis  the  PMi0  measurement  at  x=1.0  m  z=  1.39  m 
on  tower  B  (just  adjacent  to  the  road). 
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The  results  indicate  a  general  agreement  between  the  PMio  and  the  QUIC  simulations  trends,  as 
was  the  case  for  the  mean  wind  and  turbulence  profiles.  The  quiet  zones  exhibit  substantial 
removal  as  shown  by  the  bottom  two  measurements  in  each  of  the  two  top  panels.  The  PMio 
concentration  measurements  in  the  displaced  profile  (Region  C  in  Fig.  8)  are  significantly  higher 
than  those  in  the  quiet  and  bleed  zones.  This  is  a  result  of  two  mechanisms:  1)  The  windbreaks 
are  tremendously  effective  at  removing  PMio  from  the  flow  field  and  2)  the  very  low  porous 
windbreaks  permit  a  very  low  flow  rate  to  pass  through  them;  to  compensate  the  flow  contains  a 
component  that  flows  over  the  windbreaks.  This  flow  component  is  present  at  the  windward  side 
of  the  first  windbreak.  Because  this  is  also  the  location  of  the  road  PMio  is  transported  over  the 
windbreak  and  cannot  be  removed.  The  flow  of  PMio  over  the  windbreaks  can  be  greatly 
mitigated  if  a  windbreak  is  installed  upwind  of  the  road  along  with  the  downwind  windbreaks. 
The  upwind  windbreak  strategy  has  the  undesirable  consequence  that  concentrations  near  the 
road  will  be  significantly  higher  than  they  would  without  the  upstream  windbreak. 

In  an  effort  to  better  quantify  the  performance  of  the  model  with  the  limited  data  sets,  the 
following  perfonnance  metrics  were  computed  as  done  for  Handford,  WA  above:  (FB),  NMSE, 
and  RE.  Here,  Mean  Wind  is  a  comparison  of  mean  wind  speed,  Turbulence  is  a  comparison  of 
the  local  u*  and  Concentration  is  local  concentration. 


Quantity 

FB 

NMSE 

RE 

Mean  Wind 

0.12446 

1.29956 

-38.46 

Turbulence 

-0.32409 

0.55432 

-19.84 

Concentration 

-0.14765 

0.14106 

-9.62 

Table  9:  Error  metrics  for  the  Raft  River  Valley  Windbreak  experiment  for  the  mean  wind, 
turbulence  and  concentration  measurements. 


Figure  41  below  illustrates  the  model  performance  using  scatter  plots.  The  RE  for  the  wind  and 
turbulence  field  for  windbreak  (Malta,  ID)  case  is  larger  than  for  the  deep  vegetation  field 
experiment  (Handford,  WA),  however  in  spite  of  the  larger  flow  error  the  concentration  errors 
are  similar  in  both  cases  (12%  for  Handford  and  9.6%  Malta). 
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Malta  Mean  Wind 


Malta  Turbulence 


Malta  Concentration 


Figure  41:  Normalize  scatter  plots  comparing  the  QU1C  predicted  mean  wind  speed,  turbulence  and 
concentrations  with  the  field  observations  for  the  Raft  River  Windbreak  field  experiment. 
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4.3.  Development  of  Simple  Empirical  Formula  for  Native  Vegetative  Removal 

A  set  of  many  QUIC-PLUME  simulations  was  perfonned  to  estimate  canopy  removal  as  a 

function  of  the  height  and  thickness  of  the  native  vegetation.  The  simulations  were  run  with  the 

same  basic  QUIC  set-up  described  in  Section  4.2.1.  The  results  are  shown  in  Figure  42  together 

with  results  from  four  field  studies  including  the  Handford,  WA  field  study.  Fluxes  were 

computing  by  counting  particles  as  they  passed  horizontal  planes  over  a  time  interval.  These 

studies  include  Dugway,  Utah  (Veranth  et  al.  2003),  Ft.  Bliss,  Texas  (Etyemezian  et  al.  2004) 

and  Las  Cruces,  New  Mexico  (Speckart  et  al.  2013).  The  QUIC  simulation  results  and  the  four 
•  •••  •  •3^  (  # 
field  studies  indicate  increases  in  canopy  height  ( H  )  and  canopy  density  Tm  exponentially 

decrease  the  transmitted  fraction,  TF.  For  relatively  sparse  and  short  canopies,  shown  on  the 

lower  left  hand  corner,  a  relatively  small  increase  in  either  density  or  height  greatly  reduces  TF 

and  for  tall  and  dense  canopies,  on  the  upper  right  hand  corner,  TF  is  insensitive  to  changes  in 

either  canopy  height  or  density.  This  suggests  that  there  is  a  “law  of  diminishing  returns”  in  the 

mitigation  of  fugitive  dust.  Once  a  critical  height  and  density  are  reached,  the  addition  of  more 

vegetation  would  not  enhance  the  canopy’s  removal  ability. 


It  should  be  noted  that  Figure  42  also  indicates  the  sensitivity  of  the  model  to  changes  in 
parameters  related  to  H*  and  T Regions  where  the  contour  lines  are  packed  closely  together  are 
particularly  sensitive.  Recalling  Eq.  (43) 


f*  _  OtU-xHcan 

im  —  1 


mcany 


(43) 


We  see  that  T^  includes  key  parameters  related  to  leaf  area  (via  a),  roughness  (via  the  canopy 
heights),  atmospheric  stability  and  turbulence  (via  K  and  u*).  Hence,  this  non-dimensional  result 
provides  key  insight  into  the  sensitivity  of  a  wide  range  of  variables  effecting  deposition  in 
canopies. 
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Figure  42:  Contours  of  constant  transmitted  fraction,  TF,  as  a  function  of  H  and  7^  together  with  the 
field  data  from  the  Dugway  (  | ),  Ft.  Bliss  (A),  Las  Cruccs(^^)  and  Handford  (  Q  )  field  studies. 
Contours  are  generated  from  the  results  of  56  QU1C  simulations. 

Lastly,  the  56  QUIC  simulations  can  be  fit  by  the  following  equation: 

TF  =  (1  -  CF)  =  (1  -  exp(— 2.8/T))  exp(-2.0:C0-64)  +  exp(-2.8/T).  (48) 

In  Eq.  48,  CF  and  TF  are  the  captured  fractions  and  transmitted  fractions  respectively. 

This  equation  is  compared  to  the  experimental  data  in  the  plot  shown  in  Figure  43.  The  ordinate 
of  the  plot  is  physically  interpreted  as  the  transmission  efficiency  (1 -removal  efficiency)  and 
incorporates  the  effects  of  varying  the  canopy  height  for  a  given  plume  height.  Canopies  at  the 
far  left,  such  as  Ft.  Bliss,  are  sparse  and  inefficient  at  removing  particles.  Dense  canopies  such  as 
Dugway  are  efficient  and  intermediate  canopies  such  as  Las  Cruces  are  intennediate  between  the 
two  extremes.  Eventually,  as  at  a  canopy  density  similar  to  that  of  Dugway  (7^~3),  the  removal 
efficiency  reaches  a  limit  and  increases  in  canopy  density  fails  to  increase  removal  efficiency. 
The  fitted  model  to  the  56  QUIC  plume  simulations  slightly  over  predicts  removal  compared  to 
the  experimental  data. 
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TF  —  exp(—2. 8H*) 
1  —  exp(— 2. 8  //*) 


m 


Figure  43:  Eq.  34  compared  to  the  four  field  studies  measuring  ambient  removal.  This  figure  shows  that 
it-...  •  •  .  TF-exp(-2.8H*)  . 

the  PM|0  transmission  emciency,  - - - —increases  with  increasing  lm  .  Consequently  the 


removal  efficiently,  1  — 


1- exp  (-2. 8//*) 
TF—exp(—2.8H*) 


l-exp(-2.8H  ) 
Increases  with  increasing  T ^ . 
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5.  Conclusions  and  Implications  for  Future  Research/Implementation 


This  SERDP  project  produced  a  wide  range  of  important  results  that  have  improved  our 
understanding  of  flow  and  deposition  in  vegetative  canopies  and  around  windbreaks  as  well  as 
improved  modeling  capabilities.  Specifically,  we  have  (i)  shown  experimental  evidence  that 
turbulence  enhances  dry-deposition  onto  all  surfaces  of  vegetation  even  for  small  particle  sizes 
(e.g.  PMio),  (ii)  developed  a  dry-deposition  model  for  the  removal  of  particles  that  incorporated 
the  effects  of  turbulence,  which  were  neglected  in  previous  dry-deposition  models,  (iii) 
developed  accurate  models  for  the  mean  wind  and  turbulent  flow  resulting  from  windbreaks,  (iv) 
conducted  a  field  study  at  Handford,  WA  that  demonstrated  the  utility  of  real-time  QUIC 
modeling  for  experimental  design  and  prediction  of  transport  within  native  vegetation,  (v) 
developed  a  simple  empirical  relationship  and  methodology  that  relates  PMio  removal  to  native 
vegetative  height  and  thickness  (vi)  conducted  a  second  field  study  in  the  Raft  River  Valley  near 
Malta,  ID  in  which  mean  and  turbulent  flows  were  measured  along  with  PMio  transport  within 
and  downwind  of  windbreaks  and  produced  data  used  to  validate  the  new  windbreak  and  particle 
models  integrated  into  QUIC.  The  results  from  the  enhanced  QUIC  model  show  good  agreement 
with  the  field  data  from  the  Handford,  WA  and  Malta,  ID  studies  along  with  data  from  the  other 
field  studies.  Specifically  the  predicted  concentration  relative  error  at  Malta,  ID  was  -10%  and 
the  error  at  Handford,  WA  was  -12%. 


Implications:  (i)  The  improved  QUIC  model  can  be  run  to  estimate  the  effectiveness  of 
vegetation  or  windbreaks  on  near-source  removal  in  the  field,  minimizing  the  need  for  costly 
field  experiments  or  providing  a  way  to  help  design  experiments,  (ii)  The  Simple  Empirical 
Formula  that  was  developed  can  be  very  useful  for  practitioners  to  make  first  estimates  of  the 
effectiveness  of  vegetation  in  remove  near-source  dust  at  a  specific  location.  The  non- 
dimensional  fonn  of  this  method  provides  insight  into  the  sensitivity  of  a  wide  range  of  variables 
effecting  deposition  in  canopies,  (iii)  The  new  deposition  model  that  was  produced  yields 
plausible  PMio  removal  rates  while  assuming  realistic  vegetation  element  sizes.  Previous  models 
required  unrealistically  small  element  sizes  for  plausible  PMio  removal  rates.  This  implies  that 
measurements  of  typical  vegetation  length  scales  (e.g.  leaf  width)  can  be  used  in  the  new  model. 

Future  research:  This  work  has  improved  our  understanding  of  how  turbulence  influences  PMio 
deposition  to  vegetation.  However,  to  be  more  practically  useful,  more  generalization  is  needed. 
Further  field  experiment  data  are  necessary  to  better  understand  the  new  deposition  model’s 
ability  to  be  widely  applicable  over  different  types  of  vegetation  and  windbreak  designs.  For 
example,  Finch  (1988)  provides  a  list  of  species  of  vegetation  that  are  particularly  tolerant  to 
mechanical  damage  from  wind  blown  soils  that  should  be  better  understood  in  the  context  of  the 
new  models.  Further,  the  USDA  (1999)  provides  recommendations  and  standards  for 
windbreak/shelterbelt  establishment.  Understanding  the  effectiveness  of  these  recommended 
windbreaks  in  tenns  of  turbulent  transport  and  deposition  on  a  species-specific  basis  is  a  critical 
need.  In  addition,  the  QUIC  model  should  also  be  evaluated  for  a  larger  range  of  particle  sizes  so 
as  to  ensure  its  usefulness  for  other  practical  applications  including  windbreak  interactions  with 
snow  and  sand  particles.  Additional  field  datasets  would  be  necessary  to  do  this. 


67 


6.  Literature  Cited 


Arya,  S.  P.,  2001.  Introduction  to  micrometeorology,  2nd  Edition.  Academic  Press,  New  York. 

Amatul  U.,  2006.  Implementation  and  validation  of  particle  transport  model  for  vegetation  in 
Quick  Urban  Dispersion  Modeling  System.  MS  Thesis,  Department  of  Mechanical  Engineering, 
University  of  Utah,  Salt  Lake  City,  Utah. 

Aylor,  D.,  Ferrandino,  F.,  1989.  Dispersion  of  spores  released  from  an  elevated  line  source 
within  a  wheat  canopy.  Bound. -Lay.  Meteorol.,  46,  251-273. 

Aylor,  D.  E.,  Flesch  T.K.,  2001.  Estimating  spore  release  rates  using  a  Lagrangian  stochastic 
simulation  model.  J.  Appl.  Meteorol.,  40,  1 196-1208. 

Breda,  N.J.J.,  2003.  Ground-based  measurements  of  leaf  area  index:  a  review  of  methods, 
instruments  and  current  controversies.  J.  Exp.  Botany,  54,  2403-2417. 

Bradley,  E.  F.,  Mulheam,  P.  J.,  1983.  Development  of  velocity  and  shear  stress  distributions  in 
the  wake  of  a  porous  shelter  fence,  J.  Ind.  Aerod.,  15,  145-156. 

Bruun,  H.H.,  1995.  Hot-wire  anemometry:  principles  of  signal  analysis.  Oxford  University  Press, 
New  York. 

Chapra,  S.C.,  Canale,  R.P.,  2006.  Numerical  methods  for  engineers.  McGraw  Hill,  New  York. 

Cheng-I,  H.,  Katul,  G.,  Tze-wen,  C.,  2000.  An  approximate  analytical  model  for  footprint 
estimation  of  scalar  fluxes  in  thermally  stratified  atmospheric  flows.  Adv.  Water  Resour.,  23, 
765-772. 

Choudhury,  B.  J.,  and  Monteith,  J.  L.  (1988).  A  four-layer  model  for  the  heat  budget  of 
homogeneous  land  surfaces.  Q.  J.  Roy.  Meteorol.  Soc.,  1 14,  373-398. 

Cionco,  R.M.,  1965.  Mathematical  model  for  air  flow  in  vegetative  canopy. 

J.  Appl.  Meteorol.,  4,  517-522. 

Comte-Bellot,  G.,  Corrsin,  S.,  1965.  The  use  of  a  contraction  to  improve  the  isotropy  of  grid¬ 
generated  turbulence.  J.  Fluid  Mech.,  25  657-682. 

Counihan,  J.C.R.,  Hunt,  Jackson  P.S.,  1974.  Wakes  behind  two-dimensional  surface  obstacles  in 
turbulent  boundary  layers,  J.  Fluid  Mech.,  64,  529-563. 

Etyemezian  V.,  Ahonen  S.,  Nikolic  D.,  Gillies  J.,  Kuhns  H.,  Gillette  D.,  Veranth  J.,  2004. 
Deposition  and  removal  of  fugitive  dust  in  the  arid  southwestern  United  States:  measurements 
and  model  results.  J.  Air  Waste  Manag.  Assoc.,  54,  1099-1  111. 


68 


Finch,  S.J.,  1988.  Field  windbreaks:  design  criteria.  Agr.  Ecosyst.  Environ.,  22,  215-228. 

Guan,  D.,  Zhang,  Y.,  Zhu,  T.,  2003.  A  wind-tunnel  study  of  windbreak  drag.  Agr.  Forest 
meteorol.,  118,  75-84. 

Hinds,  W.C.,  1982.  Aerosol  technology:  properties,  behaviors,  and  measurement  of  airborne 
particles.  Wiley,  New  York. 

Judd  M.J.,  Raupach  M.R.,  Finnigan  J.J.,  1996.  A  wind  tunnel  study  of  turbulent  flow  around 
single  and  multiple  windbreaks,  parti:  velocity  fields.  Bound. -Lay.  Meteorol.,  80,  127-165. 

Mao,  Y.,  Wilson,  J.D.,  Kort,  J.,  2013.  Effects  of  a  shelterbelt  on  road  dust  dispersion.  Atmos. 
Environ.,  79,  590-598. 

May,  K.R.,  Clifford,  R.,  1967.  The  impaction  of  aerosol  particles  on  cylinders,  spheres,  ribbons, 
and  discs.  Ann.  Occup.  Hyg.,  10,  83-95. 

Marple,  V.,  Liu,  B.,  1974.  Characteristics  of  laminar  jet  impactors.  Environ.  Sci.  Technol.,  8, 
648-654. 

Moran,  S.M.,  Pardyjak,  E.R.,  Veranth,  J.M.,  Speckart,  S.O.,  Edmunds,  B.L.,  Stoll,  J.R.,  2011. 
Developing  strategies  for  fugitive  dust  mitigation  and  transport  flux  using  native  vegetative 
windbreaks  for  dust  control.  Proceedings  for  Student  Poster  and  Paper  Session.  Air  and  Waste 
Management  Association,  Orlando,  FL,  Paper  #16. 

Pace,  T.G.,  2005.  Methodology  to  estimate  the  transportable  fraction  (TF)  of  fugitive  dust 
emissions  for  urban  scale  air  quality  analyses.  USEPA.  (8/3/2005  Revision) 
http://www.epa.gov/ttnchiel/emch/dustffactions/transportable_ffaction_080305_rev.pdf 

Pardyjak,  E.R.,  Brown,  M.J.,  2001.  Evaluation  of  a  fast-response  urban  wind  model:  comparison 
to  single  building  wind-tunnel  data.  Proceedings  of  the  3rd  International  Symposium  on 
Environmental  Hydraulics.  D.  Boyer  and  R.  Rankin  (Eds.),  Tempe,  AZ. 

Pardyjak,  E.R.,  Speckart,  S.,  Yin,  F.,  Veranth,  J.M.,  2008.  Near  source  deposition  of  vehicle 
generated  fugitive  dust  on  vegetation  and  buildings:  model  development  and  theory.  Atmos. 
Environ.,  42,  6442-6452. 

Perera,  M.,  1981.  Shelter  behind  two-dimensional  solid  and  porous  fences,  J.  Wind  Eng.  Ind. 
Aerod.,  8,  93-104. 

Pope,  S.B.,  2000.  Turbulent  flows.  Cambridge  University  Press,  New  York. 

Raupach,  M.R.,  Woods,  N.,  Dorr,  G.,  Leys,  J.F.,  Cleugh,  H.A.,  2001.  The  entrapment  of  particles 
by  windbreaks.  Atmos.  Environ.,  35,  3373-3383. 

Seinfeld,  J.H.,  Pandis,  S.N.,  1998.  Atmospheric  chemistry  and  physics:  from  air  pollution  to 
climate  change.  Wiley,  New  York. 


69 


Singh,  B.,  Hansen,  B.,  Brown,  M.J.,  Pardyjak,  E.R.,  2008.  Evaluation  of  the  QUIC-URB  fast 
response  urban  wind  model  for  a  cubical  building  array  and  wide  building  street  canyon,  Environ. 
Fluid  Mech.,  8,281-312. 

Singh,  B.,  Pardyjak,  E.R.,  Norgren,  A.,  Willemsen,  P.,  2011.  Speeding-up  urban  fast  response 
Lagrangian  dispersion  simulations  using  video  game  technology,  Environ.  Model.  Softw.,  26, 
739-750. 

Slinn,  W.G.N.,  1982.  Predictions  for  particle  deposition  to  vegetative  canopies.  Atmos. Environ., 
16, 1785-1794. 

Speckart  S.O.,  Pardyjak  E.R.,  2013.  Removal  of  PMio  in  the  near-source  zone  downwind  of 
unpaved  roads  part  1 :  field  study,  in-prep,  to  be  submitted  to  Atmos.  Environ. 

Tabler,  R.D.,  1991.  Snow  fence  guide.  Strategic  Highway  Research  Program,  Report  SHRP- 
W/FR-9 1-106.  National  Research  Council,  Washington,  D.C. 

Tabler,  R.D.,  1992.  What's  new  in  fences  for  snow,  sand,  and  wind  control,  Proceedings  of  the 
97th  Annual  Conference,  American  Railway  Bridge  and  Building  Association,  Chicago,  Illinois, 
188-195. 

Tennekes,  H.,  Lumley,  J.L.,  1972.  A  first  course  in  turbulence.  The  MIT  Press,  Cambridge. 

T.S.I.,  2010.  Model  8520  DustTrak  Aerosol  Monitor  Operation  and  Service  Manual  1980198, 
Revision  S,  Shoreview,  MN. 

USDA,  1999.  Natural  Resources  Conservation  Service.  CORE4  Conservation  Practices.  Chapter 
3j :  Windbreak/Shelterbelt. 

Van  Ulden,  A.P.,  1978.  Simple  estimates  for  vertical  diffusion  from  sources  near  the  ground. 
Atmos.  Environ.,  12,  2125-2129. 

Venkatram,  A.,  2004.  On  estimating  emissions  through  horizontal  fluxes.  Atmos.  Environ.,  38, 
1337-1344. 

Veranth,  J.M.,  E.R.  Pardyjak,  Seshadri,  G.,  2003.  Vehicle-generated  dust  transport:  analytic 
models  and  field  study,  Atmos.  Environ.,  37,  2295-2303. 

Warner,  S.,  Platt,  N.,  Heagy,  J.F.,  2004.  Comparisons  of  transport  and  dispersion  model 
predictions  of  the  URBAN  2000  field  experiment.  J.  Appl.  Meteorol.,  43,  829-846. 

Williams,  M.D.,  Brown,  M.J.,  Boswell,  D.,  Singh,  B.,  and  E.R.  Pardyjak,  2004:  Testing  of  the 
Quic-Plume  model  with  wind-tunnel  measurements  for  a  high-rise  building,  Fifth  AMS 
Symposium  on  the  Urban  Environment,  Vancouver,  BC,  J5.3. 


70 


Wilzack,  J.M.,  Oncley,  S.P.,  Stage,  S.A.,  2001.  Sonic  anemometer  tilt  corrections  algorithms. 
Bound. -Lay.  Meteorol.,  99,  127-150. 

Woods,  N.,  Craig,  I.P.,  Dorr,  G.,  Young,  B.,  2001.  Spray  drift  of  pesticides  arising  from  aerial 
application  in  cotton.  J.  Environ.  Qual.,  30,  697-701. 

Zhu,  D.,  Kuhns,  H.D.,  Gillies,  J.A.,  Etyemezian,  V.,  Gertler,  A.W.,  Brown,  S.,  2011.  Inferring 
deposition  velocities  from  changes  in  aerosol  size  distributions  downwind  of  a  roadway.  Atmos. 
Environ.,  45,  957-966. 


71 


Appendix  A.  T.S.I.  DustTrak  Colocation  Tests 

The  DustTrak  sensors  used  in  both  of  the  new  field  experiments  reported  here  were  rented  from 
Galson  Labs  (http://galsonlabs.com/)  and  Argus-Hazco  (http://argus-hazco.com/).  The 
companies  calibrate  the  instruments  prior  to  shipment.  DustTrak  monitors  are  calibrated  to  PMi0 
concentrations  by  the  utilizing  the  respirable  fraction  of  Arizona  Road  Dust  (ISO  12103-1,  Al). 
Further,  just  prior  to  each  of  the  field  tests,  on-site  colocations  test  were  conducted.  Each  test 
consisted  of  placing  all  of  the  instruments  at  2.5  m  above  ground  level  and  sampling  for  95  min. 
All  DustTraks  were  located  2  -  3  m  downwind  of  traffic  generated  by  a  vehicle.  Table  Al  shows 
the  results  from  the  Malta  Raft  River  Windbreak  Experiment  colocation  test.  The  data  indicate  a 
standard  deviation  of  0.32  mg/m3  if  DT  RS1  is  included.  If  that  sensor  is  omitted  the  standard 
deviation  is  0.068  mg/m  .  For  the  Hanford  tests  (Table  A2)  the  standard  deviation  was  0.062 
mg/m3. 


DT_1 

0.413 

mg/m3 

DT_2 

0.432 

mg/m3 

DT_3 

0.413 

mg/m3 

DT_4 

0.384 

mg/m3 

DT_5 

0.536 

mg/m3 

DT_6 

0.482 

mg/m3 

DT_7 

0.554 

mg/m3 

DT_8 

0.436 

mg/m3 

Table  A2  -  Comparison  of  the  concentrations  measured  during  colocation  tests  for  the  Handford, 
WA  tests. 
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Appendix  B.  List  of  Scientific/Technical  Publications 

1 .  Peer-Reviewed  Journal  Articles: 

a.  Moran,  S.M.,  E.R.  Pardyjak,  and  J.M.  Veranth,  2013.  Understanding  the  role  of 
grid  turbulence  in  enhancing  PMi0  deposition:  Scaling  the  Stokes  number  with  R>, 
Physics  of  Fluids,  25,  115103;  http://dx.doi.Org/10.1063/l.4829451. 

b.  Speckart,  S.O.  and  E.R.  Pardyjak,  A  method  for  rapidly  computing  windbreak 
flow  field  variables,  submitted.  Journal  of  Wind  Engineering  and  Industrial 
Aerodynamics,  December  2013. 

c.  Speckart,  S.O.,  J.M.  Veranth,  and  E.R.  Pardyjak,  Removal  of  PMio  in  the  near¬ 
source  zone  downwind  of  unpaved  roads.  Part  I:  Las  Cruces,  NM  field  Study,  in¬ 
preparation  for  submission  to  Atmospheric  Environment,  manuscript  available 
upon  request. 

d.  Speckart,  S.O.,  J.M.  Veranth,  and  E.R.  Pardyjak,  Removal  of  PMio  in  the  near¬ 
source  zone  downwind  of  unpaved  roads.  Part  2:  Quantifying  near-source  capture, 
in-preparation  for  submission  to  Atmospheric  Environment,  manuscript  available 
upon  request. 
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a.  Moran,  S.,  2012.  The  influence  of  vegetative  windbreaks  on  fugitive  dust  (PMio) 
mitigation  and  transport  flux:  the  role  of  turbulence  in  enhancing  particulate 
deposition.  MS  Thesis,  University  of  Utah. 

b.  Speckart,  S.O.,  2013.  Numerical  modeling  and  field  measurement  of  vehicle¬ 
generated  road  dust  transport  and  removal  near  the  road.  PhD  Thesis,  University 
of  Utah. 

c.  Price,  T.,  2014  (In  preparation).  Toward  understanding  the  role  of  turbulence  in 
enhancing  particle  deposition  onto  vegetation.  MS  Thesis,  University  of  Utah 

3.  Conference  Proceedings 
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2011.  Developing  strategies  for  fugitive  dust  mitigation  and  transport  flux  using 
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